{
 "cells": [
  {
   "cell_type": "code",
   "execution_count": 1,
   "metadata": {},
   "outputs": [],
   "source": [
    "import re\n",
    "import pandas as pd\n",
    "import numpy as np\n",
    "import matplotlib.pyplot as plt\n",
    "%matplotlib inline\n",
    "\n",
    "from statsmodels.api import tsa\n",
    "from dateutil.parser import parse\n",
    "\n",
    "from sklearn.linear_model import LinearRegression\n",
    "from sklearn.metrics import mean_absolute_error\n",
    "\n",
    "def parse_quarter(string):\n",
    "    \"\"\"\n",
    "    Converts a string from the format YYYYQN in datetime object at the end of quarter N.\n",
    "    \"\"\"\n",
    "    \n",
    "    # Note: you could also just retrieve the first four elements of the string\n",
    "    # and the last one... Regex is fun but often not necessary\n",
    "    year, qn = re.search(r'^(20[0-9][0-9])(Q[1-4])$', string).group(1, 2)\n",
    "    \n",
    "    # year and qn will be strings, pd.datetime expects integers.\n",
    "    year = int(year)\n",
    "    \n",
    "    date = None\n",
    "    \n",
    "    if qn=='Q1':\n",
    "        date = pd.datetime(year, 3, 31)\n",
    "    elif qn=='Q2':\n",
    "        date = pd.datetime(year, 6, 30)\n",
    "    elif qn=='Q3':\n",
    "        date = pd.datetime(year, 9, 20)\n",
    "    else:\n",
    "        date = pd.datetime(year, 12, 31)\n",
    "        \n",
    "    return date\n",
    "\n",
    "\n",
    "alcohol_consumption = pd.read_csv('data/NZAlcoholConsumption.csv', \n",
    "                                  parse_dates=['DATE'], \n",
    "                                  date_parser=parse_quarter,\n",
    "                                  index_col='DATE')\n",
    "alcohol_consumption.sort_index(inplace=True)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "# Modelling and Forecasting Time Series\n",
    "## Getting a stationary time series\n",
    "\n",
    "In the previous notebook, you had seen that differencing with a lag of 4 made sense. \n",
    "\n",
    "Let's do that, does the resulting series look stationary?"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 2,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "image/png": "iVBORw0KGgoAAAANSUhEUgAAAe8AAAF1CAYAAADBdGLoAAAABHNCSVQICAgIfAhkiAAAAAlwSFlzAAALEgAACxIB0t1+/AAAADl0RVh0U29mdHdhcmUAbWF0cGxvdGxpYiB2ZXJzaW9uIDMuMC4xLCBodHRwOi8vbWF0cGxvdGxpYi5vcmcvDW2N/gAAIABJREFUeJzsvXl4XPV59/29Z180o32XbMub8AoGQ1lCQoBgaAi42ReatUnavOn7Pk8SN9CnydU3aQsJbdKW8LRJmo20Kdl4XBpIHHYCGIKNwcY2smXJkiUv0midff09f5w5o9HozMw5M0czR9L9uS5dlzRz5sxvjmbm+7t3EkKAYRiGYZilg6naC2AYhmEYRhss3gzDMAyzxGDxZhiGYZglBos3wzAMwywxWLwZhmEYZonB4s0wDMMwSwwWb4apMEQkiGh9nvs+RES/LeGcTxPRn+iwttNEdGO559FhHf9KRF+q9joYxqiweDPLkrQIhYnIT0TTRPQCEf0pEVX0Pa9VVIUQ/yGEuGkx17TYENEPiehvNBz/USJ6Lvs2IcSfCiG+qv/q9IOI6ojoASI6S0STRHQfEVG118WsDFi8meXMO4QQHgCrAdwD4IsAvlfdJTFGhYgsGh/SAOAQgM3pn1sBvEfvdTGMEizezLJHCDEjhHgYwPsAfISItgIAEdWmLadxIhoior/KtsyJ6ONEdJyIpohoHxGtTt9ORPRNIhojohkiOiyfMxsi+lsA1wL4FhEFiOhbWXffSEQn0+e+X7bYsq1Qtc+j8Lym9GsZSj/2ASKqzbr/NiI6mvZIPE1Em/Kc5yIiGiSi9yvcp7g2IvoUgA8B+Iv0a/7v9PF3EtGptCfkGBH9Ufr2TQD+FcBV6eOn07fPs96J6JNE1J+2cB8moo6s+0Taq6J0PdcR0ZNENEFEPiL6DyKqy3rsaSL6IhEdBhAkoj1E9Muc13ofEf1j7jUQQgwIIb4phJgWQpwHcAJAa7H/D8PoghCCf/hn2f0AOA3gRoXbhwH8Wfr3BwD8FwAPgDWQvnw/kb5vN4B+AJsAWAD8FYAX0vftAnAQQB0ASh/TnmcdTwP4k5zbBIBfpR+/CsA4gJvT930UwHPlPA+Aj6fXvhZADYCHAPw4fd9GAEEAbwNgBfAX6WNt2dcNwKXpa3VrnufLuzYAPwTwNznHvwdABySD4X3pNbTnvuas4zPnAHA9AF96TXYA9wF4VuX1XJ9+rXYAzQCeBfCPOe+TVwF0A3ACaE+vrS59vwXAGIDLirzf3g1gEkB3td/7/LMyftjyZlYaZwE0EJEZkojcJYTwCyFOA/gHAH+cPu7TAO4WQhwXQiQA/B2AS9LWdxyS4F8EgNLHnNO4jnuEZLENA3gKwCUKx5T6PB8C8A0hWYYBAHcBeH/aLfw+AI8IIR4TQsQB/D0k0bo66/HXAngYwEeEEL/K8xya1iaE+LkQ4qwQIiWE+CmAkwCuUPFa5NfzfSHEK0KIaPr1XEVEa7KOUbyeQoj+9GuNCiHGAXwDwFtyzv/PQogzQohw+jU8izn3980AfEKIg/kWR0RvAvBdALcJIc6ofE0MUxYs3sxKoxOShdQEwAZgKOu+ofT9gBQn/6e0a3k6/RgC0CmEeBLAtwDcD+ACEX2HiLwa13E+6/cQJAt5HmU8T4fC67JAcunOu08IkQJwBnOvGwD+FJKX4al8T6B1bUT0YSJ6Net6boX0P1BD7poDACZy1qx4PYmohYgeJKJRIpoF8O8Kz5sruD8CcEf69zsA/LjI+j4D4JtCiOeKHMcwusHizawYiOhySF/4z0Fyw8YhibTMKgCj6d/PAPi0EKIu68cphHgBAIQQ/yyEuAzAFkiu6D15nrassX0aniebs1j4uhIALuTel44Nd2PudQOSeK8iom+WuLZ5rzntrfgugM8CaBRC1AF4HdJmaMHxxV4PEbkBNOasOR93p8+/XQjhhSTGuRnhuc+/F8D2dH7BrQD+o8hztKfXyDAVg8WbWfYQkZeIbgXwIIB/F0IcEUIkAfwMwN8SkSctMJ+DZJkBUhLVXUS0JX2OWiJ6T/r3y4noD4jICik+GgGQzPP0FyDFnktZt5bnyeY/AfxPIuohohpILv+fpt3/PwPwdiK6IX3ezwOIAngh6/F+SO7iNxPRPSWsLfc1uyEJ5Hj6sR+DZHkj6/guIrLleT0/AfAxIrqEiOzp1/NSOtRRDA+AAIBpIuqEis2PECIC4Bfp5/192hVfiHehuMAzjK6weDPLmf8mIj8kK/p/QYp3fizr/j+HJDwDkKzxnwD4PgAIIf4PgK8BeDDtbn0dwC3px3khWZJTkNy5E5Bix0r8E4B3p7Og/1nj+rU8Tzbfh+TqfRbAICRh/fP06+qDZH3eB8n78A5IJXWx7BMIIaYhJXrdQkRK9daF1vY9AJvTLvK9QohjkPIJ9kMS6m0Ans8615MAjgI4T0S+3CcSQjwB4EsAfgngHIB1ABZkwOfh/4eU6DYD4BFIyXtq+FF6ncVc5oAk3O9VeV6G0QUSoiyvHsMwzLKDiFYBeANAmxBittrrYZhc2PJmGIbJgqRa/88BeJCFmzEqWjsKMQzDLFvSyXAXIIUCbq7ychgmL+w2ZxiGYZglBrvNGYZhGGaJweLNMAzDMEsMw8a8m5qaxJo1a6q9DIZhGIapGAcPHvQJIZqLHWdY8V6zZg0OHDhQ7WUwDMMwTMUgoqHiR7HbnGEYhmGWHCzeDMMwDLPEYPFmGIZhmCUGizfDMAzDLDFYvBmGYRhmicHizTAMwzBLDBZvhmEYhllisHgzDMMwzBKDxZthGIZhlhiG7bDGMAxjJPYeGsW9+/pwdjqMjjon9uzqxe4dndVeFrNCYfFmGIYpwt5Do7jroSMIx5MAgNHpMO566AgAsIAzVUEXtzkR3UxEfUTUT0R3Kty/ioieIqJDRHSYiP5Qj+dlGIapBPfu68sIt0w4nsS9+/qqtCJmpVO2eBORGcD9AG4BsBnAB4hoc85hfwXgZ0KIHQDeD+B/l/u8DMMwleLsdFjT7Qyz2OhheV8BoF8IMSCEiAF4EMDtOccIAN7077UAzurwvAzDMBWho86p6XaGWWz0EO9OAGey/h5J35bNXwO4g4hGADwK4M+VTkREnyKiA0R0YHx8XIelMQzDlM+eXb1wWs3zbnNazdizq7dKK2JWOnqINyncJnL+/gCAHwohugD8IYAfE9GC5xZCfEcIsVMIsbO5uegscoZhmIqwe0cn7n7nNlD6266zzom737mNk9WYqqGHeI8A6M76uwsL3eKfAPAzABBC7AfgANCkw3MzDMNUhN07OmExEZxWM56/83oWbqaq6CHeLwPYQEQ9RGSDlJD2cM4xwwBuAAAi2gRJvNkvzjDMkiGWSCGeFAjHk4gnU9VeDrPCKVu8hRAJAJ8FsA/AcUhZ5UeJ6CtEdFv6sM8D+CQRvQbgPwF8VAiR61pnGIYxLKFYIvN7IJIocCTDLD66NGkRQjwKKREt+7YvZ/1+DMA1ejwXwzBMNQjG5uq8/ZEE6t22Kq6GWelwb3OGYRgVhLMs79lIvIorYRgWb4ZhGFUEo/Mtb4apJizeDMMwKghmWd5+tryZKsPizTAMo4IQW96MgWDxZhiGUQFb3oyRYPFmGIZRQSjGljdjHFi8GYZhVBCMZlneURZvprqweDMMw6hAtrwb3DbMhtltzlQXFm+GYRgVhGJJ2Cwm1Lms7DZnqg6LN8MwjApCsQRcNjM8Dis3aWGqDos3wzCMCoLRJNw2C7wOC1veTNVh8WYYhlHBnOVt4VIxpuqweDMMw6ggGEvCZbfAY+eYN1N9WLwZhmFUEIom4M5Y3izeTHVh8WYYhlFBMJaEy2aBx2FFOJ5EPJmq9pKYFQyLN8MwjApCsQTcdjO8TgsAIMDWN1NFWLwZhmFUEMqyvAFukcpUFxZvhmEYFYSic9nmALjWm6kqLN4MwzBFSKUEQvFkJmENYPFmqguLN8MwTBEiiSSEAFx2C7zsNmcMAIs3wzBMEYJRaShJtuXN4s1UExZvhmGYIoRiklDPT1hjtzlTPVi8GYZhipCxvO1seTPGgMWbYRimCNmWt9VsgsNqYsubqSos3gzDMEUIxeYsbwDwOLi/OVNdWLwZhmGKIFveTqvkMuexoEy1YfFmGIYpQnbMG5Asb67zZqoJizfDMEwRsmPeAOBxWDDLljdTRVi8GYZhihDMiXl7HVZOWGOqCos3wzBMEULRBIgAh0V2m3PMm6kuLN4MwzBFCMaScFnNMJkIgCzebHkz1YPFm2EYpgihWAIuuyXzt8dhRSSeQjyZquKqmJUMizfDMEwRQjFpopgMd1ljqg2LN8MwTBGC0SSctvmWN8D9zZnqweLNMAxThFAswZY3YyhYvBmGYYoQjCXnxbzlmd7cqIWpFrqINxHdTER9RNRPRHfmOea9RHSMiI4S0U/0eF6GYZhKEIqy5c0YC0vxQwpDRGYA9wN4G4ARAC8T0cNCiGNZx2wAcBeAa4QQU0TUUu7zMgzDVIpQLJnprgbMWd4s3ky10MPyvgJAvxBiQAgRA/AggNtzjvkkgPuFEFMAIIQY0+F5GYZhKkIwlsh0VwPmLO/ZMLvNmeqgh3h3AjiT9fdI+rZsNgLYSETPE9GLRHSz0omI6FNEdICIDoyPj+uwNIZhmPIJRedb3jXsNmeqjB7iTQq3iZy/LQA2ALgOwAcA/BsR1S14kBDfEULsFELsbG5u1mFpDMMw5RFPphBLpubFvK1mE5xWM5eKMVVDD/EeAdCd9XcXgLMKx/yXECIuhBgE0AdJzBmGYQxNKD2UxJkl3gD3N2eqix7i/TKADUTUQ0Q2AO8H8HDOMXsBvBUAiKgJkht9QIfnZhiGWVTkcaBu+/z8Xo/DAn+ULW+mOpQt3kKIBIDPAtgH4DiAnwkhjhLRV4jotvRh+wBMENExAE8B2COEmCj3uRmGYRabYFSyvF0LLG8rW95M1Si7VAwAhBCPAng057YvZ/0uAHwu/cMwDLNkyFjetoWW9yyLN1MluMMawzBMATKWt32+5e11WjlhjakaLN4MwzAFyGd5ezlhjakiLN4MwzAFCKazzd12pZg3W95MdWDxZhiGKUA4bXm7cmPedgsi8RRiiVQ1lsWscFi8GYZhCpA/21zussbWN1N5WLwZhmEKEMpnefNwEqaKsHgzDMMUIBhLwmom2Czzvy55LChTTVi8GYZhChCKJhZY3UC25c1uc6bysHgzDMMUIBhLzhtKIpMZC8qWN1MFWLwZhmEKEIol4LIvtLy9bHkzVYTFm2EYpgDBaGHLm2PeTDVg8WYYhilAKJYv5s3izVQPFm+GYZgChGLJBTXeAGAxm+CymdltzlQFFm+GYZgChGJJxZg3IE8WY/FmKg+LN8MwTAGC0YRizBvgmd5M9WDxZhiGKYDkNs9vebN4M9WAxZthGCYPQggEY4kFE8VkeLIYUy1YvBmGYfIQiacgxMK+5jJseTPVgsWbYRgmD8H0UJJ8lrfXYeEOa0xVYPFmGIbJQygzDjSf5c1uc6Y6sHgzDMPkIRRPW975ss3tFkQTKcQSqUoui2FYvBmGYfIRTFvezjzi7XVyf3OmOrB4MwzD5CGUiXnnT1gDuEUqU3lYvBmGYfIQzMS885eKASzeTOVh8WYYhslDxvIuUCoGgFukMhWHxZthGCYPwVja8s7bpEV2m7N4M5WFxZthGCYPoWhhy9ubdptzrTdTaVi8GYZh8iBb3k5rMcubxZupLCzeDMMweQjHEnDZzDCZSPH+Gju7zZnqwOLNMAyTh2AsmTfTHAAsZhNcNjNb3kzFYfFmGIbJQyiayNsaVUYaTsKWN1NZWLwZhmHyUMzyBqSkNba8mUrD4s0wDJOHUCyRt7uaDI8FZaoBizfDMEwegtHiljdPFmOqAYs3wzBMHkKxRN4abxkPz/RmqgCLN8MwTB6C0WTe7moybHkz1UAX8Saim4moj4j6iejOAse9m4gEEe3U43kZhmEWk3A8WdTy9rLlzVSBssWbiMwA7gdwC4DNAD5ARJsVjvMA+H8BvFTuczIMw1SCYDShIuZtQSyRQjSRrNCqGEYfy/sKAP1CiAEhRAzAgwBuVzjuqwC+DiCiw3MyDMMsKolkCtFESkWdN48FZSqPHuLdCeBM1t8j6dsyENEOAN1CiF8VOhERfYqIDhDRgfHxcR2WxjAMUxqhuGRJu4vGvLm/+d5Do7jmnifRc+cjuOaeJ7H30Gi1l7Ts0UO8lZr+isydRCYA3wTw+WInEkJ8RwixUwixs7m5WYelMQzDlEYomh4HqtryXplJa3sPjeKuh45gdDoMAWB0Ooy7HjrCAr7I6CHeIwC6s/7uAnA2628PgK0Aniai0wCuBPAwJ60xDGNkgrH0ONAilrd3hVve9+7rQzg+P94fjidx776+Kq1oZaCHeL8MYAMR9RCRDcD7ATws3ymEmBFCNAkh1ggh1gB4EcBtQogDOjw3wzDMosCWtzrOToc13c7oQ+F3pQqEEAki+iyAfQDMAL4vhDhKRF8BcEAI8XDhM6ws9h4axb37+nB2OoyOOif27OrF7h2dxR/IMExFyVjeKrLNAazYcrGOOidGFYS6o85ZhdWsHMoWbwAQQjwK4NGc276c59jr9HjOpYgcG5JdTHJsCAALOMMYjHAsbXkX6W3uTVves+GVaXnv2dWLL/z8NSRSmVQnOK1m7NnVW8VVLX+4w1oF4dgQwywdZMu7WJ13zQqPee/e0YntnbWZvzvrnLj7ndvYIFlkdLG8GXVwbIhhlg5zMe/C4m02Edw284oVbwCIJlMAgFu2tuFf7risyqtZGbDlXUHyxYA4NsQwxmMu5l3cxlnJ/c1TKYFT4wEAK9f7UA1YvCvInl29cFrn7+I5NsQwxiSUiXkXtryBlT3T++xMGJG4ZHnPrtANTDVgt3kF2b2jE0II/M+fvQYAaKqx4a/evpljQwxjQILRBCwmgs1c3MbxOCzwR40lXJWqbOkfk6zuVq99xW5gqgFb3hXmhs2tmd+5TIxhjEsoloTLZgaRUhPJ+Uhuc+MIVyW7nsnivaO7fsVm3FcDFu8K4/NHM7+PTHGiGsMYlWA0AXeRMjEZr9NY4l3JypZT40HUuaxY3eSCP5KAEKL4g5iyYfGuMBPBWOZ3Fm+GMS6heLJoprmMFPM2jtVZycqWU2MBrG+ugddhRSw9iY1ZfFi8K4xsedc6rTgzGaryahiGyUcomijaGlXG47AYqsNaJStbTo0HsK65Bl5nulmNgTYxyxkW7wrjC0jifXF3HVveDGNggjH1lrfXYUUskUIkx1VdLSpV2TIVjGEiGMP6lprMgJbZsHE2McsZFu8KMx6Q3OYXd9Xigj+CaMIYH3aGYeYTiqmPeRttpvfuHZ24+53bMvOaW732Rel6Jtd3r2txZ9rEGil8sJxh8a4wE4Eo6l1WrG50Qwjg3HSk2ktiGEaBUFRbzBswlnDdur0dcurY/R+8dFHLxNY3e+B1ruwBLZWGxbvC+AJRNNXY0V0vxZ7Ydc4wxiQYS6jqrgYAHrtsdRpHuLJF9Pzs4hgJp8YDsFlM6Kx3rvjRqJWGxbvC+AIxNNXY0dXgAgCcmeKkNYYxIqFoUlV3NcB4bnMAmMmquT4/szji3T8WwNomN8wmypquZpxrsJxh8a4wE4Eomjx2tHrssJgIIyzeDGM4hBDaLG8DWp3ZDVMuLJLl3T8ewPqWGgDZc82Ncw2WMyzeFcYXiKHRbYPFbEJ7nYPd5gxjQKKJFFJCXV9zYAlY3rPRAkeWRiSexMhUGOuaJfF22cwwm8hQG5jlDIt3BYnEkwhEE2j22AEAXXUuFm+GMSCZoSRW9aVigLGsTlm8G902XFgEt/nAeBBCIGN5E5FU785u84rA4l1BxtMNWppqbACA7gYnN2phGAMSjEoC5FJZKlZjYMt7Y6tnURLWMmViacsbkDYxbHlXBhbvCiK3Rm10py3vehfG/FHDNHZgGEZCtrzVxrzNJkKN3VhjQWXx7m2TxFvvnuP9YwEQAWub3ZnbjNZpbjnD4l1B5NaoTbLbPF0uthj9hhmGKZ1gTLa81bnNAeP1N58Nx2GzmLCqwYVYIoXpkL5rOzUeQFe9E46s0AJb3pWDxbuCyK1RZbd5V71ULsZxb4YxFqGoNssbkK1O4wjXbCSOWqcVbbUOAPrXevenB5JkwzHvysHiXUFkt3lTjWR5dzdIljfXejOMschY3io7rAHGm+k9E5bEu9Wrv3gnUwKDvmAmWU1GGo1qnA3McobFu4KM+6OosVsybqYWjwNWM7HlzTAGI5QWb7W9zQHZbW488c5Y3jpmnI9OhRFNpOYlqwGS25xj3pWBxbuCSK1RbZm/zSZCR52TxZthDMZcwppWy9s4VudMOA6vw4IWjx1E+oq3nGmea3l7HBYEogkkU/omxzELYfGuIBPp1qjZdNU7ucsawxgMOebt1CTexrS8rWYTGt12XbusyQNJFlje6ZneAQNdh+UKi3cF8QWiaMyyvAGgu96FM5NseTOMkZiLeS9ht3lIEm8AaKu16xrz7h8LoNFtQ717/vcZt0itHCzeFUSeKJZNV70TvgDXejOMkQjFknBYTTCbqPjBabwOK2LJlCE+y6mUgD+amBNvr0N3t3mu1Q0Ys9PccoXFu0LEkylMheIK4s3lYgxjNIJR9UNJZLwG6rLmjyYgxJwbu9Xr0M1tLoRA/3gA61qUxDtteXO52KLD4l0hpuQyMc9CyxsAx70ZxkCEYurHgcoYabKYPFEs2/KeCsV18QpMBmOYDsWxLquzmoy8WTDCNVjuaNtaMiUzLjdoyYkRlWJ57z00inv39eHsdBgddU7s2dWL3Ts69Vsss2Lh95ZEKZa3kSaLya1RM5Z3ulxsbDaKVY2uss4tJ6vlZpoD2THv6l+D5Q6Ld4XwBZQt7xaPHTazSXWjlr2HRnHXQ0cQTu+gR6fDuOuhIwCwIr9kGf3g99YcoVhSU4MWYM7yNkK8d0bB8gakRi3livep8SAAZfH2Gsj7sNxht3mFmMi0Rp0v3iYTobNefa33vfv6Ml+uMuF4Evfu69NnocyKhd9bcwRjCU0NWgBjWt5z2eb6dVnrHwvAaTWjo9a54L4ajnlXDBbvCiH3Nc8tFQPkWm914p1viAkPN2HKhd9bc4RjSThVzvKWmRPv6luduTFvuUWqHnO9T40HsLbZDZNCJr7VbILLZjbENVjusHhXCF8gBpvFBI/Cbr6r3olRlW7zjrqFu91CtzOMWvi9NUdplrfsMq6+1ZlreXsdFjitZt0sb6UyMRmjDWhZrrB4VwhfIIrmGjuIFu5Wu+pd8AVimX7KhdizqxdW8/xz2Mwm7NnVq9tamZXJnl29MOe8P+2WlfneCkW1x7xr7MZJ1poJx2ExUeY1EBHaah1li3colsDodFgx3i3jNdiAluUKi3eF8AVi8/qaZyOXi42qcJ3v3tGJHd11MBFAQOYDeuPmVj2Xy6xArl7XiJQQcNvNkCX84q7aFZesBpRmeZtNhBq7MWZ6z4Tj8Dqt84yFVq+9bLf5QDpZrZDl7XVa2fKuALqINxHdTER9RNRPRHcq3P85IjpGRIeJ6AkiWq3H8y4lfP4oGnOS1WS0lovNhBO4rrcFg/e8HT/99JWYicRx72/e0G2tzMrkJ78fhgDwqz+/FoP3vB0fvXoNXhmexmgJMe+9h0ZxzT1PoufOR3DNPU9i76FR/Re8SCRTApF4SrPlDRRvkVqp6yL3Nc+mzVu+5Z1vIEk2i9Umdim/pxaDssWbiMwA7gdwC4DNAD5ARJtzDjsEYKcQYjuAXwD4ernPu9TInSiWTbeGRi2hWAInx/zY1lkLALhsdQM+ctUaPPDiEF4+PanfgpkVRSyRwn+8NIzrepvR0yQ13/jkm9cCAL777ICmc8klZ6PTYQjMlZwtlS/bzDhQjXXegOwyVrY6K3ldZMs7m9ZaB8ZmoxCi9Ilfp8YCMBGwpil/uZnXYc0kzOnFUn9PLQZ6WN5XAOgXQgwIIWIAHgRwe/YBQoinhBCyMr0IoEuH560Y5e74UimByeDCiWIyzR477BYTzqiwvI+enUVKANu7ajO37dnVi45aJ774y8OG6KvMLD1+c/Q8xv1RfOSqNZnbOuuc2L2jEw++PJwpdVTDUi85k8eBau2wBhS2Oit5XWYjCUXLO5ZMYTLd7bEU+scDWNXggt2S/9oshuW91N9Ti4Ee4t0J4EzW3yPp2/LxCQC/VrqDiD5FRAeI6MD4+LgOSysfPXZ8M+E4EimR121OJNd6F7e8D4/MAEDG8gYAt92Cu9+5DQPjQdz35EnV62IYmQdeOI3VjS68ZWPzvNv/9C3rEE2k8IPnT6s+11IvOQtGS7e8CwlXJa/LbB63OVBerfepsWBBlzkwF/Mux8LPRe9rtxxc8HqIt9LYHcX/GhHdAWAngHuV7hdCfEcIsVMIsbO5uVnpkIqjx47Pl2nQouw2B6S4t5qY95GRabR5HWhJfxBl3ryxGe+6tAvffmYAR8/OqF4bw7w+OoMDQ1P44ytXL6jdXd9Sg5u3tOFH+0+rTsSqc1kVb18qJWey5a1llreMp4DbvJKleFLMe/7mQ26RWuqAkkQyhUFfsGCyGiBtYOJJKW9AL/S8dsvFBa+HeI8A6M76uwvA2dyDiOhGAP8LwG1CCPU+uCqjx45Pbo3anMfyBtQ3ajk8OoNtWS7zbL506ybUuaz44i8PI5HU74OzlFgOO+pK88D+03BazXjPzm7F+z9z3Xr4Iwn8+4vDRc914oIfgUgCuf07nFbzkik5k8W7VMs7X6nYbRd3LLhtMa6LEEKKeTvyWN4zpX39jkyFEUumFKeJZbMYLVL37OqF3TJfrkq9dsvFBa+HeL8MYAMR9RCRDcD7ATycfQAR7QDwbUjCPabDc1YMPXZ8Gcvbk1+8u+tdmAzGMi47JfyROAbGg9jeqSzedS4bvnL7Vrw+Oot/e25Q9fqWC8tlR11JpoIx/NerZ/FHl3YucLPKbOuqxbUbmvC95wYL5lQEogn86b8fRK3Lhr9+x+ZMtnZnnRN3v3Pbkik5C6YT1kqLeUuWd67L+MxkCP/58jDavHa0py1gs4kW5boEY0lQh8DrAAAgAElEQVQkU2LB/7PZYwdR6W5zeSCJGssb0LfH++4dnXhH1uanrdZR8rVb6mEdmbLFWwiRAPBZAPsAHAfwMyHEUSL6ChHdlj7sXgA1AH5ORK8S0cN5Tmc49uzqXdAmUeuOL9Ma1V3IbS5nnOd/A70+OgsAeS1vALhlaxt2bWnFNx87gYF0WcdKodI76uVg5f/0wBlEEyl8+KrC1ZufuW49fIEofn7gjOL9Qgh88ReHMTQRwrc+uAMfvroHn37zOgDAk194y5IRbkBq0AKUbnnHkwLRxJznKxRL4JMPHEAqJfDgp67C/rtuwJ5dvUimBK7f1KLbumVyu6vJWM0mNNWUXuudKRMrIt5ylrvezWpSWRui//zklSW/p5ZLJ0FdpooJIR4F8GjObV/O+v1GPZ6nGshvkC/+8jCiiRSsZu275YlADGYTod6lRrxD6G3zKB5zZHQaALC9qy7veYgIX719K274xjP4kx8dQCSRxLnpyIoY71jJHXWlJ3AtxqjOZErgx/uH8Ac9DbiozVvw2CvXNmDHqjp8+9kBfOCKVbCY5+/7v//8aTxy5BzuvOUiXLm2EcDce/rsdCRTfrYUyFjeJcS8vVlWp8NqhhACf/GLw+i74McPPno51qSvw+Z26Xq/cc6PK3oadFq5xExIWbwBoL2MLmv9YwE01dhRmyenQSZzDXQuF3t1eBoOqwmReKosl/yeXb2485eHEcnaYJkJ+MJNG/VYZsXgDmsq2L2jM1PmZTaRYuyqEL5AFA1um2Ijfxk1jVoOj8ygq96JhgIWPAC0eB14+7Z2DPiCODsdWTEu5OY8YYnF2FHns/L/5pFjiKfzDfSyzBcrHPDE8QsYnQ7jo1evKXosEeEz163HyFQY/314fkrLgdOTuPvR43jb5lZ8Ol0bDsyJ95lJdX37jUJIzjbX2GENWNjf/DvPDuBXh89hz65eXNc7Z2Vv7pDE+9giJJfK7mol8W71OkpOWOsfD2B9S/FNmDfnGujBVDCGAV8QV6U3huVMLdu9oxOfv2nOc+p1WJAUUthnKcHirYJIPImzM2F01jkRiadUz96W8QWiBV3mgJSJ7rCaCn7RHRmdmVffXYjfnVxYamf0pIxyxC6eTMFmXvh2XqxEqXzWvC8Qw6VffQy3fet32POL13QR3MUKBzywfwjttQ68TWVr3RsuakFvqwf/8vQppFKSC9MXiOL/+ckr6Kx34u/fc/G8dpzdDdo6BxqFoFznXWKHNUASrmdPjONrv3kDb9/Wjj97y7p5x7V47Gh023Ds3Gz5C85BdpvnNmkBSu+yJoTAqSIDSWQWY675q2ckr6O8ASr33LJ36Lsf3olXv3wT3trbjK/+6jiOjCydSh0WbxUMT4YgBDJfcicuaIsljwdiea1CGSIqWC42E4pjaCKEbZ35XebZnJ1W/oAaNSmjXOvy28+cwsh0GB+9enUmq9brsCxKQtB0SAqDKNHgtuLt29rx+ugs4sn5SUulCu5ihAP6x/x4rt+HO65cvcAFng+TifBn163DiQsB7PzbxzObLJ8/iv/9oUsXWHqtXgcsJlLVv8BIhGIJmE20ILtZDa+NSCKz+/7n8ZEf/B6tHju+/u7tCwYSERE2d3hx9OziibeS5d1W68B0KK65mdN4IIrZSKJojTcAeJ36zzU/NDwFEwHXrG8CUL5LXna7exwWmEyEb7z3EjTV2PCZnxzMXD+jw+KtgkGf1Iz/pox4+zU9fiIQzdtdLZuueidGppW/6I6MSjtCtZb3UkvKKMe6fOP8LP7piZN4+/Z2/PVtW/HiX96AS1fVobvBpbtw+yNxfOQHL0MIscDSd1rN+PKtW3DPu7YjX38KrYI7HYrBYlbeKJTzv3xg/xBsZhPed7lyeVg+UqkUCMBkMAYBIJpIwWQinFTY0JpNhI469bPqjUIoloTLalacAFiIvYdG8e1n5lrJCgFMhuJ47NgFxeM3d3hx8kIAsYS+ZZ2zBSzvzFxvjdb3qTHpO1CNeDutZlhMpGvM+9CZaVzU5kWrV/oeLdfy9qdd5PIkuHq3Dfd98FKcm47gL37xmq4NZhYLFm8VyOK9pbMWnXVOTeIthCjY1zybQrXeh9PJals71Im3HlnylaRU6zKeTOELP38NXocVX7ltS+b2W7a24+jZWQxP6Gf1hWNJfOKHB/D66Az+5Y7L8PV3b0dnnROEheVQ+YS1WL5CNtOhGO743ktIpZQ3CqX+L/2ROH55cAS3bm9XtanM5h8eO7mgA1M8KfJusrpUdg40EqFosqQysXv39c3LMgekzU2+a7O53YtYMpXJ4taLmXAcRIBHIWY/V+utTbz7x9WViQGSV0HPFqmplMCrw9PYsaoObpsFJirfqg+kH59dC3/Z6np88eaLsO/oBU0dBasFi7cKTvuCaHTbUOu0YmNrDfrOqxfvYCyJSDyVtzVqNl31LkyH4oqZlIfPzGBNo6topqfM7h2duPud2zI71cVyIetFPrErFm7416dP4fXRWfzN7q3zrvHNW9sAAL9+/Zwu64vEk/jUjw/gwNAk/vF9l+CmLW3YvaMTz995PQbveTuev/P6eddWafNEACaCMXz7mVNFd/YzoTju+N5LOHE+gH/76OX4+ru3Z6yEjjJqXAHglwdHEIwl8REViWq5aN1kqW0+ZCSCsURJZWJar82WTNKavq5zuUGLUoJsW630GdEa9z41FoDLZs7UqBdDz7Gg/eMB+KMJXLqqHiYTwaPD4BP5O7bGMf///CfX9uBtm1vxd48ex6HhqbKeY7Fh8VbBoC+YKfHY2OrBwHhQdQeziUxr1OLi3V0g4/zI6Ay2FSgRU2L3jk689Jc3YlO7Fxe1ew0r3IAkdkph5JlwDE/3Kff1OX5uFv/85Encur0dt2xrn3dfd4ML2zpr8evXz5e9tngyhc/+5BX87qQPX3vX9nnNIvIhb56yLfOvvWsb3r6tHXf/+g187mev5Y07Zgv3t//4Mry1twW7d3Ti3ndvBwD8yx2Xlfy/TKUEHtg/hIu763Bxt7b3E6A9HNNV78KYP7qkBuaEYqVZ3lqvTU9TDRxWk+5Ja0p9zWVKdpuPS8lqakMJelresojuWCW9X73O/F3s1BLIcZvLEBH+/t0Xo63WgY/94GVcdfcThu3jwOKtgtMTQaxpnBPvWDKF0yrdsWr6msvka9TiC0QxOh3O21mtGDduasHBoSlMlTFNaLF560UtEEL6MMli9+VbN2Ftswcf/+HLC6xV2V1e67TiK7dvVTznLdva8OqZ6ZISu7Iz37f/9T48fnwMX719S94WokrkWubvvXwVvvXBHfj82zbi/xwaxXu/vX+B+3ImHMcff/8l9J3341//+FK89aK58iK5Oc/hUe0ZsfLrWfuXj2LAF8TWTuVeAsXQGo6R39OlzASvFsFoAq4SLG+t18ZsIlzU5tV9FoE0DlR5/R6HFW6bGec0us1PjQVUxbtl9BwLemh4GnUua6ZXgB7n9kcSsFtMsCkkJda6rHjvzi5Mh+M4N2PcUltdmrQsZ0KxBC7MRtGTnl+7sVX60jt5wa/qzTzulwRTbcIasHCut5ysVqizWiFu3NSK+57sx9MnxvBHO4w5jfXpvjEIAD/6+OW4bPVc04r3X7EKe35+GHf/+g28cd6Pq9c14h8fP5kRg49dszpvHPmWre34+m/68JvXz+Pjb+pRvZbcBizheAqWtLuuXIgIf37DBmxs8+BzP30V7/jWc/jQH6zCzw+M4Ox0GBYzIZUS+O5HduL6i+aXcHXWSTX+R0amARTuiFbo9QDALw+O4vLVjZotePl4tQ1jssvF1MRLjUAollS12c5F67UBpKS1X712FkIIzQly+ZgpYHkD0oASLZZ3MJrA2ZkI1jWrb7TjcVhw2qdPrsMrw1PY0V2XuT5S//jyE9Y8jvzy99OXRxbcJifQGsWDyeJdBPkNKLvN17fUgAjou+Bf4KpVwqfBbd7gtsFpNS+wvI+MzIBoLkamlW2dtWjx2PH4MeOK92PHLqCpxo5Luuvn3e6yWfCtD+7Apqc8+PvfnsDeV0fnZXI/+PsRXNxVr/iB6mly46I2D379+jlN4q2U+Z5ICV0/uLu2tOGhz1yDD3x3P/7x8bkxrvGklJym1ISCiLCtszYzFlYtSq8nEk+V/Hp27+hU/bh8G1IjE4wlsMruKumxWq4NICWt/eSlYYxOhzONmsplJhxHW4HYdJvXoSlhbWBcfaa5jNehT8x7NhLHybEAbt0+F6ryOqwYKjMRNRBJLHCZZ7MU+p+z27wIpyekN67ssnHazFjV4FIsjVFiIj1RrFHFTp6I0N3gXNCo5fDIDNY2uUu2/Ewmwg2bWvDMiXHNZSmV6N8dTSTxdN843ra5RbF+mojw2es3oMFtW1CCVayc7Jat7TgwNIUxDZZGpT64vW0e2C0LY6uxZP4M5e1dtTg5FkA4pj6GXM0vohaPA1YzLamktVA0CXcJDVpKYfMiJK3NhBMFLe82rwMXZtVPFusflxJ0tYi3HkllgJSoK8RcvBuQkuHKnVjmj8QLfp8uhVJbFu8iyGVicswbkFznfSrLxXyBKOpcVlhVNsJQatRyZHQaF2tMVsvlxk2tCEQTeGlwQvVjKjWla/+pCQSiCdy0ua3gcfli9oVE6A+3tUEIYN9R9Ylrlfzg5rOA8r2mbZ21SKYEjp1Tb31X84toKdZ6h2KlxbxL4aI2D4igW9KaEAKz4bhijbeM7DaXu+QVYu+hUXxp71EAwIe/93vVn32v04JgLFn2aOJDw1MgwrzkSsmqLz9hrZDlvRRKbVm8izDoC6LFY5/X57i31YNBXxDRRHHrR01r1Gxy62IvzEZwYTZacrxb5pr1TXBYTXg8T8MIJSo1peu3xy7AbTPjqnWNBY8rRYQ2tHqwrtmtKev8CzdtRK79v1gfXK2vSR5Ko8V1Ln0R6VcnrpWlVOsthJCyzStkebtsFqxtcuvWaS2aSCGWTBW1vBMpgYkiCazy5l3OzD47E1G9eZet2nL7hb8yPIUNLTXz6rG9TgsC0URZGwN/JLGgTCwbuVqkNp341+YtrzxzMWDxLsLprDIxmQ2tNUimRMYqL8REIKapEUZXvROzkUSmRZ/ca1dtZ7V8OKxmXLuhGY8fH1PdPagS7tZUSuCxYxfwlt5mOKyFvzBL3Q3/4bZ2vDgwkSnbK0azxwEBoM5lVWzAoidaX1NbrQMtHrumHsy7d3Tif9y4IfN3pedrd9Xlb/trNGLJFBIpUdJQklLZ3FGrm9s809e8gEtYbblYOZt3r6P8FqlCCBw6M40dOXkwXh02Bv5I4YQ1QPrc/OsdOwEAX3v3dkMJN8DiXZTTE0H0NM4Xb3lkp5pmLb5AFE1FGo1kM1frLVkqh0dnYCJgc3t54g1IJWOj02EcP6fO5V8Jd+urI9MY90eLuswB5dppNSJ089Y2pIRk4avhe88NoKnGhhfvukGxAYuelPKatnfVai4Xq3VK3p8nP/+WRX09SnTVOzG+RGq95VnelbK8ASlpbXQ6nBnlWQ6F+prLyMlsxZLWytm8y5Z3OX3CT0+EMB2Kz4t3S+eWR46WLt6BaEKxA10u8ujW44swQKZcONu8AP5IHL5AbIHl3dPkhjlPP+dcxgNRXKvJbT5XWrOloxZHRqaxsdUDpw5fJtdf1AqiI3ji+IVMokwhvnDTRnzu56/NSxLT29362LELsJgIb80al1gIrdm8gPQBXN3owq9fP48PXLGq4LH9YwE81TeO/3HjhqKeAL3Q+pq2ddbhiTfGisbtsjk4NIX6rFrZSpJdLqYl6akayLO8S+mwViqZpLVzs0VDR8VQJd5yi9QilndHnVOxPl/N5l2P4SSvDEnNWS5dnWN5O8ubWiaEkMRbRQJwrcuKjlqHIcWbLe8CyGVico23jN1iRk+Tu2jSWjSRhD+S0Ow2B6QvOiEEDo/MYFuJzVlyafbYcUl3HR4/rs4CbfU6IAQycR/PIrRY/e3R87hybaPqtq+lQES4eWsbXuj3FbVufvD8IGwWE+64Un0ddaXZ3lULIYCjGqzvg8NTuGx1vW61xFpYSuViIXkcaAkd1kpFtu70SFqT39+FxLupxgYTFXeb79nVC0tO9Yfazbvs2i6nXOzQmSl47Basz+kPUO65w/EkkilRMOadzeYOr+4tbPWAxbsAg5kysYXWwsbWGpwsIt5ymZgWt3mdS+qANDIVwtmZCCaCsbLj3dncuKkVr43MqGrS8E9PnESr146X/vJGdNU7cf1FLboKd/9YAKfGg7hpi7p50uXwh1vbkUgJPFZg4zIVjOGXr4xg9yUdmgd2VJKt6c3cEZXiPRmMYWA8uMCCqRRdBdr+Go1gtPKWd7PHjhaPXZdOa2osb4vZhGaPvajbfPeOTnTWSaV+WnM/ZIEtx/I+NDyNi7vrFvRol636Ut3m8lAStV6rTe1eDPiChgv7sHgX4HQ6IW1148LmCRtbPRiaDBWst5UbtGjJNpfnep+ZDKc7aUFzT/NC3LhJEson31DuFy7z0sAEXhqcxKffvA4Oq+RpUJOgpwV5VKK8psVke5c0Ee7XR/IPKvnJ74cRiac0NXSpBs0eOzpqHaozzuXe0Jetqo54t3jsS6bWO2N5VzDmDehn3cnWaCHxBtKNWops4APRBEamI/jktWs1537MCWxp1nEolsAb5/0L4t1A+Za3XGZWLGFNZlO7F8mUUN3bo1KweBdg0BdER61DMfa5sdUDIVBwnF+mu5oGyxsAuhuk0prDIzOwmAgXtZXWh1qJja016G5wFi0Zu+/JfjTV2DMx4rVNbgyOB3Wdc/vbY+exrbO2IvXGsuv8dyd9ig0eYokUHth/Gm9a34SL2krrZFdJtnXV4nB6c1eMg0NTsJgoU2ZWaUwmQmfd0igXk8W7ktnmgOQ67x8LqCo/LYRseRcTplZv8RapL5+eRDIlSorDy1ZtqZb34ZEZJFMClypsODPiXeLGQM5S1yLegPGS1li8CzCoUCYmI/c4L5Rx7ku7zZs1umC76l0YnQrj8MgMets8uiZOERFu3NSK5/p9CMWUP1gHhybxXL8Pn37z2kyiXE+TG/5oIvOaymVsNoJDw9O4afPiW90yt2xtQyyZUvQ6PHrkHC7MRvEJg1vdMtu76nB6IqQqQ/ng0BS2dHh1SXosFaXmQ0ZE/kxU+lpt7vAioYN1NxOOo8ZugaVIU6i22uItUl88NQGrmbAza9aAWixmE9w2c8nW8aFhaWN6icLku5oyy9Dm3Obq8mxWN7jgspl1n/5WLizeBTg9kV+81zS6YDObcGKskHin3eYahxx01TvhjyZwcGhK13i3zI2bWhFNpPDcSZ/i/f/8RD8a3DZ86Mq5zOyedNKI3C62XB4/LgnoTVuKl4jpxaWr6tHisePXR+Y3bBFC4HvPDWJtsxtv2dhcsfWUg5zE+HqROGk8mcJrI9NVi3fLLJVGLcF0qVglY94AsKVD+n+WKxDFhpLItHodmI0kCob99g9MYEd3fckbGY+j9Damh4an0NPkRr1CyNFsInjspQ8nkdek1vI2pb2fLN5LhOlQDNOh+IIabxmL2YS1zW6cKGR5+2Nw2cyaWy3K2bnheHJRXJ1X9DTA47DgieMLLdDXzkzjmRPj+JNre+ate216EzM4ro94//bYeaxudGFja+VKh0wmyXX+9ImxeV6Hl09P4cjoDD5+Tc+C5BijIot3sbj38XOziMRTuKzK4t3d4IIvENPUk70ayO+LSmabA1nWXZlx72KtUWWKlYvNhON4fXQGV5ZRuuZ1WkpKKhNC4JXhacV499y5rSUnrPnzzPIuxKZ2L46fm9U1bFguLN55yPQ0L1AXu7HVgxMF3FwTwWhJWcsnx+bO+Y3HTujeS9xqNuG63hY88caFBf2N73vyJOpcVnz4qjXzbu+oc8JmNmFAh6Q1fySOF/oncNPm1oqXLt28tQ2ReApP941nbvvecwOoc1nxrkuNOXFNiXq3DasaXDgyWjjufTBdK1tt8Z6b621s61u2vF0VqvGXMZkIm9rLT1qTLO/iolSsUcvvByeREsBVa0sXb0+Jk8VGpsLwBaLYUSDBspyxoAGNCWuAJN7+SMJQc+lZvPOQO01Mid42D0anw3ldQ75AVPNc4L2HRnH/U/2Zv8f90UUZBnLjphb4AjG8mpX09ProDB4/PoZPXNOzYFdqNhFWNbow6Cs/4/KZE+OIJVMVdZnLXLGmAY1uGx5NZ50PT4Tw22MX8MErVlU1JlwKUtJaYcv74NAUOmodaK+t7jQkWbzPGDzuHYolYLeYisaMF4PN7V4cOzeramBIPmaLTBSTKdYidf+pCdgspoLWbzG8DktJcelDZ6TvpB0K8e65c5c+tcyvsVQMyE5aU9edshKweOdh0BeCiYBVDfln7G5Id4vKtpSz8fljaNRoed+7rw+R+PyG+4sxDOS6jdL4zSey6p6/9WQ/PA4LPnLNGsXH6FUu9tujF9Dotilmki42FrMJN21pxVNvjCEST+IHLwzCTLTA07AU2N5Zi5GpMCYLDJh4ZWgKl63RnnCkN0ul1jsYS1Q801xmc4dXKs8q4xrNhOMF+5rLZCzvfOI9MIHLVtWXlSxbquV9aHgKTqu5YJWN12kpebJYIBqH02rWtEHLTH8zULMWFu88DPqC6Kx3wmbJf4nkHuf5mrWU4jav1OzlWpcVV6xpwOPHpLj3G+dn8Zuj5/Gxa3ryfvjXNrlxeiKEZBmWQSyRwlNvjOGGTcqzuytBrdOKYCyJi770G/zw+dO4uKs282W2lJAnzeVr1nJ2OoyzMxFcVob1pBfNNXbYzCbDJ62FopWbKJbLlkyb1NKbtahNWKuxW1Bjtyi6zaeCMRw/N4ury2zV6nWWZnm/MjyN7V21BcXVW0YynNQaVdsGzW23YE2j21DlYizeeTjtC86b4a1Ed70LDqsJfecXWt7JlMBkMIZmjW7zSs5evnFzK/ou+DE8EcJ9T/ajxm7Bx/NY3YBkeccSqbI2Ei8NTsCvYnb3YrH30Ch++MLpzN8CwOtnZ3UPS1SCTKe1PPXec/Hu6lveJhOhs974c71DsWTFM81lNrZ6YDZRyeNBY4kUwvGkKvEGgFavXdFt/tLgBACU3Wfdk3Ztq03y2ntoFFff/QReOzONo0U+k1LCWulNWtS2Rs1mU7sHx8+zeBsaIQRO+4JFhziYTIQNLR6cVCgXmwzGkBLaG7RUcgi8/KF6871P4ZHD53Dl2gbUufJvNuTrUY7r/LdHL8BpNeNNG5pKPkc5KIUloomU7mGJSuB1WLG2yY3X8sS9Dw6l3Y/t+jX5KYeueidGJo1teQdjiarlPjisZqxrdpfsms20RlU5J6CtVrnL2gunJuC0msuudPE6rEikxILPmxLy7PCzaU9AIJoomOvjdVjgjyZKyg8IRNRNFMtlU5sXQxOhsmeU6wWLtwITwRj80URRyxuQZnsrNWqZa42qTbxLHXuplb2HRvEPvz0x77bn+n0Fd7s9zaWL995Do7j6nifw4xeHICDwm9fPF3/QIlCpsESl2NZVm3e29yvDU7i4uxbWKiRfKbEUGrWEYkm4K1wmlo2ctFYKavqaZ9PqdeCCgtt8/6kJ7FxTXzBkqIZMi1QV7m2ts8O9TiuEAAJ5Gk0Vwh+Jq5oolouctNZnEOvbGJ9qgyH3NJfFqhC9rR6M+aOYDs1PGsoMJdHoNgckAX/+zusXdZa00oclEi9sgTbX2FFjt2gW78yuejqSeZ7FyKBXQyXDEpVgW2ctzs9GMJZjQYViCRw9O1v1ErFsuuqdmAjG8nb2MwLBaEJzXwY92dzhxbmZSMEkxHzIIqmmzhsA2msdGPNH51mv4/4oTo4FcPW68j1jnsxwkuLirXVTPTfTW7vrXMso3Ww2yTkJBklaY/FWQBanfA1astmYTlrLrfcuta95pSjFAiUi9DS5Ndd6a91VLyaVDEtUAtm1mZu0JveGNpp4A8Coga1vKeZdPctb7rRWSmKUbHmryTYHpEYtiZSALxjN3PbigD7xbmkdlvS6im/WtG6qy5laFigx5t1R60Ct04pjBikXY/FW4PREEBYTZb5sCiH3OD+Rk3GeEW+NbvNKUaoFKpWLaav1NpKrulJhiUqxpcMLEy3stCYnq+3oNpJ4G79cLBRLwFWlUjFgzjVbynjQ2RLc5gBwYWZOvPcPTKDGbsHWjvKH82ixvLVuqmXvQimWtz+iPdsckIyXTe0ew2Scs3grMOgLorvBpaoOsKPWgRq7RUG8Y7CZTZm4j9Eo1QLtaXJjZCqsafqR0VzVlQhLVAq33YL1LTULLO9Xhqawrlm5N3S16E5vho1cLhaMVtfyfvbEOEwE/N2jb+Cae57UFFrSGvNWqvXef2oCV/Q06NKkpjYT8y5uHe/e0Ymv3L4l83exTfXcWFBtlncqJRCIlZawBkibq77z/rLKZfWCxVuBQV8IaxRmeCtBRNjQWqNoeTfW2Cre/lMtpVqga5vdEELqTKaWPbt6YTXPvw5L2VVtNLZ11uHwyEymekAIgYPDU4ZymQNAU40dNovJsJZ3KiUQjierFvOWc0NkXRidDmvKDZEnzKkW75z+5udnIhj0BctqiZqNFssbmJsg9k/vv6ToprrUeeGheBJCoCS3OSCJdziexJBOA5rKQRfxJqKbiaiPiPqJ6E6F++1E9NP0/S8R0Ro9nncxEEJgqMA0MSV6Wz3oO++fV88otUY1pstcphQLVM7A1xL33r2jE5vaPDARloWr2mhs76qFLxDNfAkP+IKYDsVLGuW4mJjSoagzBrW85byMamWbl5sbMhOWOoepzRJvrLHDbKJMxvn+AWnKoB7xbiB77rY663g4XUZYqKuljCdjeWsT77mJYtqzzQGpGgAof/qbHpQt3kRkBnA/gFsAbAbwASLanHPYJwBMCSHWA/gmgK+V+7yLxZg/ilAsWbTGO5sNrR5MhTadcAoAACAASURBVOLzZl3LlvdyQ97UnNYg3vFkCoO+EN53efeycFUbDbnTmhz3luPd1R4DqoSRy8WCmVne1bG8y80NmY2o664mYzYRWjz2zKZv/6kJ1DqtGYEqF4fVBIuJVAusNvEubaZ3oIS+5tmsb6mB2USGiHvrYXlfAaBfCDEghIgBeBDA7TnH3A7gR+nffwHgBjKoPzmTaa7R8gbmt0mdCMQMb3mXQq3TiqYam6ZysYNDU/BHE3jLxpZFXNnKZXO7F2YT4XC609rB01Ooc1kzY1yNRJeBu6yFMrO8q2N5l5sbMhOOa86xafU6Ml3W9g9M4A96GnQbi0tE8DrVtzEdmgjBbTOjQUWehtVsgstm1uw2ny1holg2DqsZ65trDDGgRA/x7gRwJuvvkfRtiscIIRIAZgDo45vRGdmiVNOgRUaeSd2XFm8hxLIVbwCay8We6huD1Uy4Zr0h/+VLHofVjI2tnjnLe3gKl66qN+Rs8q56JyaDMQQN0qUqG9nyrlbMu9wyRrV9zbNp8zpwbiaCM5MhnJkM6+Yyl/E41M/0PjMZwqpGt+o8IW8Jg0/k7milijcAw2Sc6yHeSlc6NxVPzTEgok8R0QEiOjA+Pq7wkMVncCIIm9mkKRO62WNHncuaqfWeDScQS6ZKatCyFNA6XeyZvnHsXN1QcpyJKc72zlocGZ3BdCiG/rGA4ZLVZORyMSPNRZYJxaob885OIgUAu8WkKTdkRuU40GzaaqUua/t1rO/ORssAkaHJEFY1qP/e9TrVbwxk5tzmpX8XbWqXGulMldBIR0/0EO8RAN1Zf3cBOJvvGCKyAKgFMJl7IiHEd4QQO4UQO5ubm3VYmnYGx4NY1ejSNPGKiLCxxZPJOJebHixfy7sG4/6oqg/luZkw3jjvx3W91fl/rhS2ddViOhTHw69JH71qjFtVQ5eBy8Vkb0A1O6zJSaS3X9KBphq7ptyQ2XBcdXc1mVavA/5oAk8cl8b0bmzRtw++x6FudGcqJXBmMoTVGjyepYwcnUtYK8fylmd7V9f61kO8XwawgYh6iMgG4P0AHs455mEAH0n//m4ATwq1o2YqzOmJ4tPElNjYJpWLCSHg8y938ZaT1op/AT/TJ3lQruvlePdisj2dtPbDF07DbCJc3F1b5RUpMyfebHkXYmOrB6PTYU1jL0tym9dK31FPvTGOK9c26h5qUWt5j/mjiCZS6FaRrDZ3bkvJbvNSS8WAOfGudsZ52eKdjmF/FsA+AMcB/EwIcZSIvkJEt6UP+x6ARiLqB/A5AAvKyYxAKiUwNBFCT5P6N5DMxlYP/JEEzs9GMlnnTZ7l6TZf2yyXixXvtPZ03zjaax2ZvABmcTiRHo4zMB6EiaTpbUakucYOu8WEMwacLiZb3tUaCZqN3Lnx5Ji6boaJZAqBqHa3udxlLZZM4UqdXeaAete2nGm+Wot4O62as83l48v5Hzd77Gj22KuetKbLu1QI8SiAR3Nu+3LW7xEA79HjuRaTc7MRRBMp9DRpF5q5NqkBTARLmyi2VFjV4AJR8eli8WQKz/X78I6LOwzbrGY5sPfQKL70X0czf8eTAnc9dAQADFeSR0SGzTiXa6xdVeywJiNvdk9e8KsKgciipLavuczrWZ357nviJDx2i67vGY9Ky1tueqKmTEzG69A+09sfkYaSaAmLKrGp3bss3ObLhkymeYmWNyBZQD5/FCaCqpKHpYjDakZnnbOoeB84PYVANMHx7kXGSINf1GDUWu9gVBbv6lve3fUuOKwm9J1XZ3lrbY0KSJu+bzw2NxZ4zB/Vfdqf12FFMJZEIll4pveZyRBMBHSqmCeRObdTiqdricAGovGSa7yz2dTuwckxP2KJ4rPKFwsW7yxKqfGWaXDb0FRjx4kLfowHYmhw28re3RkZNRnnT5+QS8TKHy/I5MdIg1/UIFnexnObh2IJEEnNRaqNyUTY0CIJhBpKEe979/UhEp8vPnpv+uTEsECR0sChyRA66pyaZs97HFYkUyKTq6CGQLS0oSS5bG73Ip4UODWubUiTnlT/XWogTvuCcFhNaPU4Snr8xnSPc18gumxd5jI9TW4MjgcL7nrlEjE9drpMfow2+KUYXfUuTIXiRb/QK400lMRimBCP0syEfGTE26VevCux6Zub/lX4fz08GdLkMgeyh5Ood537SxwHmosRMs5ZvLMY9EmZ5qVmXG5s9eDkWADj/uiyTVaT6Wlywx9NzGsJmw2XiFWOpTaj3KhzvUOxhCHi3TK9rR5cmI1mBo4UohTLuxKbPtnKLSawUpmYRvF2am+RKse8y2Vtkxs2i4nF2wjsPTSKZ06M443zfs2j+GQ2tnoQiiVx/Nzssi0Tk5FDC/lc53KJ2Fsv4hKxxWapzSg3aq13MJaE20BeokwejQrXuSyOWsS7Eps+NdZxIG0EaCkTm3duDUlr/khcc1KfEhazCb2tnqpmnBvnnVpF9h4axZ0PHUYiPYtPHsUHaMvW7W2TMkSjidSyF++16Yz8QV8AV/QsnF71VN8YOmod2NDCJWKVYPeOTsOKdS7yl7TRysVCUWNZ3hvSGecnLvhx+ZrCE+Jky1uLMMnvl3v39eHsdBgddU7s2dWrc7a5PLozv3V8JlMmpi3XKOOS1+A2D0T1sbwBwGUz4YVTPvTc+ciiXLtisHijcOKGln/G+qzuRMtxolg2nfVOWM2EQYVGLbFECs/3T3CJGKNIo9sGh9V4c72DsYQharxlOuuccNvMmRr+QsyE47CZTZqT7RZ70yd7AgqViw1NqJ8mlo1XxcYgl4BOMe+9h0bxyvD0gtnrQOXKM9ltDv0SN556YwxyuPw7zw7oWnJhNMwmwupGNwYVGrUcHOISMSY/Uq238crFwrEkXAboriZDRNjQ6snMTCiE3BrVaJvluZh3cct7lcaYt9aZ3smUQDCW1CXb/N59fYgn5yfrVro8k8Ub+iRu7D00irseOpLZiU2H4rrXTBqNfOViXCLGFKOr3omRaWO5zYOxpKHc5oBUwaKmXExqjWocr4GM7KIuaHlPBlHrtGruDjfnklcn3pnWqDq4zY1QnsniDeALN21ccJvWxI2l1ihDD9Y2uXF6IoRkav4O9Jm+cVy+hkvEmPwYscuaFPM21nt2Y6sHvkAME4FoweNK6WteCSxmE9w2c0HX9vBkWHOmOSA1i7JbTKqzzfUYBypjhPJMFm8AWzqlIQ51TmvJ2bpG2IlVmp4mN2KJ1LzXyCVijBpmQnFMh+LoufORkqs79CYYS8JtOMt7ru1yIWZLGAdaKaQe5Pmt4+GJoOZM8+xzq3Wbz00UK/86GaE801jbzCrx7AmprOmR/+/azCxdrXTUORVnFBu1UYYeZJeLyR++p3mKGFOEvYdG8Zuj5wEAAtVJ9lEiFEvAZTBvUW+bPKDEX3DW9kw4nhkYZDS8BUZ3JlMCI1Nh3LKtvcRzq5/pPTfLu/z/cSUy9YthrHdqlXiu34e1Te6ShRuQdmJ3PXRknuvcyI0y9KCneU6837xRsrSf5hIxpgiFkn2qJd6xRArxpDCc5d3iscPrsKCvSMa5Ud3mgOSmzufaPjsdRiIlNE0Tm39uDZa3jm5zoPrlmSvebR5NJPHSwCSu3VBectVSa5ShB801dtTYLZmkNblE7C29LYbLemWMgxFDTKGY9MVutJg3EUmdGwu4zVMpgdmIccW7kGs7k2lejttcZcKavIHQS7yrzfJ4FWVwcGgK4XgSb9pQfoy22juxSkNEWNPkwkBavOUSsbdyvJspgBFDTMH0cAu3gUrFZDa2efDI4XMQQihuiv3RBITQ1l2tkngcFpwaV7a8h0ssE5PxOiwYUdnsZ85tbszrpJUVb3k/d9IHs4lw5drCHYwYZXqaajK13nKJ2NVcIsYUwAjJPrmEosa0vAFgY0sNZsJxjPuVM85ly9NrUPEuNHd7aDIEq5nQXlvaxq20hDXj/Y9LgcW734dLV9XpkoG4EulpcmNkKoxoIomn3+ASMaY4cojJbpG+fowQYpLHShqtzhsonnFeylCSSiLHvJUmEA5PhtBV7yp5fLK0MVA30zsQlUa+GvF/XAorWryngjEcGZ3Bm9azm7dU1ja5IQTw4sAk+i5wiRijjt07OvHOS7vQ4Lbh+Tuvr3q4KWjQmDcguc0BoC/PeNDZEvqaVxKv04pESizogwEAwxOhksvEAGljEEumEE2kih4rTxRbLvk4K1q8nz/lgxDAm8pMVlvJyOViP3x+EADwVi4RY1TSVe/EZDCWSRarJqGocWPeTTV2NLhtOJlHvJeC5Q0o9yAfngyVnGkOZM8LL+4690cSht3glMKKFu/fnfDB47Dg4q7aai9lybImLd5P9Y2js86J9VwixqhEHg1qhEZGRra8AWBDSw1OFBNvlzGFSRbM3EYtM6E4ZsLxkjPNpXMX750uE4jGl1VIb8WKtxACz/X7cPW6RljMK/YylE32MJapUAz/9erZ6i6IWTLIfRWM0CY1ZOBsc0Bq1nLyQkAxtrtkLO8c8S430xzQNhY0ENVnophRWLGqNegLYnQ6rEuJ2EoldxhLKJZc9sNYGP3orDeOeAcNnG0OABtaPfBHEzg3E1lw30w4DrOJDNdgRmZOYOdbx0OTUolpeZa3Nrf5csk0B1aweP/upA8A8GaOd5fMShzGwuhHi8cBi4kUa74rjZGzzQGpXAyAoutcbtBi1EQsb57pX7LlXU7CmjxJTZXbPJ2wtlxY0eLd3eDE6kZj9gNeChixUxazdDCbSGrYUmXLe++hUXzn2QEAwHX3Pm1Iz9FcudhC8Z4JJzICaUTmYt7zBXZ4IoSmGltZgurRYHnPRhLLqiR4RYp3PJnCiwMTuJZd5mVhhLF4zNKmM0+3tUohh37kcZHykBSjCXi924Zmj12x1tvIfc2B/HHp4cnyysSALLe5qph3nN3mS51Xz0wjEE3gWu4EVhZG7JTFLC0666treS+l0E9vq0exXGwmHDdsdzUAsFtMsJppoeVdZpkYADisyufOJZ5MIRJPsdt8qfO7kz6YCLh6HYt3OazEYSyMvnTWOXHBH0FMRZONxWAphX42tNbgxIUAUqn5GeezBre8iWhBi9RYIoWz0+GyktXynVuJwDIbSgKs0MEkvzs5ju1ddYati1xKrLRhLIy+dNY7IQRwbiZc8fyTZErAbbdkXObZGDH0s7HVg3A8idHp8Dx3s9Hd5sDCsaBnp8P/t727D46rOu84/n2kXWnX0i6SbSzbKww4GKfkBUwcIM1rkzCkTChOZ9q0k06g0wx/MO20f+CWNOlffaNxkpl2ks6ESTtxkmY6yZSYMCE1L0nbkAEFGmMoGCOwwbYkv0nIkqzX1Z7+sXdXK2lXL967u/de/T4zGq2v1mfvPZb32XPuc85DzsE2H/7N8/ubLz3yLvwba+QdYhcmZjh8crjqEqAiUr3CRi31njofmZzhj/Y/y9hUdtG+2kG99VMuac05F/iRNywuIPJmlaVAS6USsWVH3lErBwprMHg//fogOQcf0P1ukYbr7si/eZ+q4zT1sXNj7Pn6L3iq9zx/s+edfOV3rg/FrZ8dXfnlYqV7nI9Pz5LNuUDf84bFI+/CMrErq9igpSCdWL6y2FxFsWD302pE52PICv289xxtLc3s2tbZ6FMRWfM2X5bArLYbtRw41Me+g0fpH55gfVsLY1MzrGuJ8d3P3cwt2zcABDJYL5ROxNlyWYLekozzoO+uVpBOxDk7MnfeJwYv0hpr4vL21urbTsY4PbJ485pSUZw2j86VrNBTr53nlu0baImtuUkHkcBpiTXRlUrUbNq8sBSskFE+eHEaA/78E9cUA3eYXNuV4ujpuZF3WIJ3uZH3FevX0XSJpUBLrShhrRC8NW0eTicGx3lzcFxVxEQCpLszSd/weE3aLrcUzAH/+tQbNXm9Wru2q53Xz40x62WchyV4L5zafnOw+mVixbaTy0+bj+ied7j9/LVzANqcRSRAMp2126glTEvBVmJHV4qpbK54zzgswTuViDM+PcvMbA7nHCd92KCl2HZrjMmZHFPZxfXCC4pLxVqD3U+rsaaC91O959lyWYK3Xa4tUUWCItORZGB4sjia9FPUdgHc6WWcF6bOR0ISvNPeHuRjk1kGL05zcXrWl2S1fNvlt18tNTo5Q6zJSMSjE/KicyXLmM05fvHaeT64Y2NgN/AXWYsynUmyOceZZZKOLsXe23aGZinYSlzjFSgp7LRWGHmnA55FnSrZxvSEj8vEYO6DwVLBu1AONErv/VUFbzNbb2aPm1mv931RCreZ3WBmT5vZS2b2gpl9uprXvFQvnBpmZDKrEqAiAVOo612LqfM7rt9KMt5EIt4U+KVgK9HWGqO7M8mrZ/OZ2yMTM5gF/15uoXDK6GSWkz4uE8u3vXxxkqhVFIPqs83vB550zj1gZvd7f/6LBc8ZBz7rnOs1s63A/5rZQefccJWvvWIHDvXxVw//HwB/9+gRcjkX2v+8IlHT3Zl/E+97a4L3XuVv2z3HBhmbmuWfP3Mjt79ri7+NN8jOrhSvnp4beadaY75kbddSsTjJxAxvDuaDd+Hf3be2l0hai1pFMah+2vxOYL/3eD+wZ+ETnHOvOud6vcf9wFmgbsPfwlKRwpTK6QuTgawaJLJW1XLk/cgL/bS1NPMbOzf53naj7OhKcez8GDOzufzWqCHY5rkwMzAymeXE0Dib0wkScX9qp8+NvJeaNs9/yImSaoN3l3NuAMD7vuT/EDO7CWgBXq/ydVcsTFWDRNaiZEszG9paOPWWv8vFprM5Hn3xNLde10WyxZ9AEQQ7N7czM+t4c/BiKPY1h/mlO08Mjvt2vxtKPxgsMW3u3fOOkmWvxsyeADaX+dEXVvNCZrYF+A5wl3OubAkhM7sHuAdg27Ztq2m+oqgtFRGJokxn0vdd1p567RwXJmb4rRu2+tpuo+3YVMg4H2NkMhuq4D3qjbzf7+P21KVT8pWMTmZ52+VrLHg75z5e6WdmdsbMtjjnBrzgfLbC89LAj4EvOueeWeK1HgQeBNi9e7cv60a2dpRfQxrWpSIiUdTdmeSV04trVVfjkcMDXJaM84FropWkes2mdposX6DkwsQMm1LtjT6lZRVGvedGpzg9MulbshpAW0szTbZMtnkEE9aqnTb/EXCX9/gu4OGFTzCzFuCHwLedcz+o8vVWbe9tO0kuuLcS5qUiIlGU6UjSPzyBc/6s9Z6cmeWxl07zm+/cHLmtkBPxZq7c0Ebv2dHQTJs3NxntrTFeHhgB/FsmBl5N72V2WRuN4LR5tb/VDwC3mlkvcKv3Z8xst5l903vO7wIfAu42s+e9rxuqfN0V27Mrw9//9rtCUTVIZK3KdCSZnMkxeHHal/Z++spZLk7Pcsf10ZoyL9ixqZ2jp8MTvCG/XOylvgsAbPNx5J1vu/L+5lPZWaazucCvhV+tqj6KOOcGgY+VOf4c8Dnv8XeB71bzOtXasyujYC0SYBlv2dCptybY6EOlqUcO97OxvTWUxUdW4tquFE8cOUPOEfhyoAWpRJx+b3MZP0fe+bZjxf3LFypsjappcxERnxWXi/mQtDY6OcNPXznLJ9+9ZdHualFx7eYUhd1kQzPy9nZCa/NWF/ja9hIj79EIFiUBBW8RCYBMZ2Gtd/XLxZ44coapbI47ro/GpizlXNs1l6QWppE3wBXr1/m+TWk6GauYsBbFWt6g4C0iAXBZMk4qEfNl5P3I4QEyHUl2XbFot+bIuHpjW3FWITQjb2/k62em+VzblRPWCkFdCWsiIjWQ6ah+rfdbF6f5n1fP8cl3bwn8lqHV+MmLpylc3X3fPxyKHSMLI2+/73eDV9O74rR5OIq3rJaCt4gEQrcPdb3/86XTZHMuslnmMLflc9a76X1ubCrwWz4fONTHD73z+/5zp3w/11QixsXpWbKzi/f/0rS5iEgNdXeuq3ra/JHD/Wzf2MY7tqZ9OqvgCduWz4UPG4UgemFixvcPG6U7uC1UDN6aNhcR8V+mI8noVLZYo3q1zo5M8vSxQT55/dZI1W1eKGxbPtfjw0Yhaa9c8Fa2uYhIDRUzzi9x9P3jFwdwDu54d3SzzKHy1s5B3fK5Hh820ksUJxmdzNLS3ERrLDrFaUDBW0QCorDW+1Kriz1yuJ+3b06xoyvl52kFTti2fK7Hh42lipOMTc1EbsocFLxFJCDm1nqvfkR2cmicX50YjlwFsXLCtuVzPT5slJYcXWh0Mhu5KXOocntUERG/bGhrIRFvWtW0+YFDfew7eLQY8FsjVoSkkjBt+Vw4z30Hj9I/PMHWjiR7b9vp6/kXa3pPlElYi2BFMVDwFpGAMDMyFUr4llPIYi5NhvrywVfZ0NYamsC2VtT6w0Zx2nwNjbzXxsdUEQmFTOe6FW/UErYlU1I7qdYYZpQtTjI6laW9NVobtICCt4gEyGpG3mFbMiW10+TVC6+UsKaRt4hIDXV3Jhm6OM34dPkiE6XCtmRKaqvS/uaaNhcRqbFuL+N8JaPnvbftJBGf/xYW5CVTUlupRGxRwppzLrIJawreIhIYc2u9lw/ee3ZluPcjb5v3d4O8ZEpqK51cPPKeyubI5lwk13lH74pEJLQKa71XmrSWjOffwnr+8mN0pRM1Oy8JvnQivihfohDMUxGrKAYaeYtIgGxKJYg12YqT1nqOD3HVhnUK3EI6uThhbaywr7mmzUVEaqe5ydjakVzRRi25nOPZN4a4+eoNdTgzCbpyCWuFoiS65y0iUmMrXS72yulRLkzMcPP29XU4Kwm6dCLG2FSWnFfnHObKgSrbXESkxjKdyRUVJ+k5PgjAzds18pZ8wppz+U1ZCoojbwVvEZHaynQkOTs6xXQ2t+Tzeo4N0d2ZLGaoy9pWLE5Sct971JtGTythTUSktro7kzgHAxcqT5075/jlG0PcdLWmzCUvncyPrkdLtkgtTJvrnreISI0VS4MukbTWe3aMoYvT3KJkNfGUKws6pmlzEZH66O5YB8CpJZLWeo4V7ndr5C15xcpipdPmU1kS8SbizdELddG7IhEJtc2XJTBbeqOWnuNDbE4n2LZ+XR3PTIKsWNN7cn7CWhQrioGCt4gETEusia5UouK0uXOOnuND3Lx9PWZW57OToCqXsDY2Fc2iJKDgLSIB1N2ZpG+4/HKx4+cvcm50SpuzyDyFID06b+QdzXKgoOAtIgGU6ay8UUvP8SFA97tlvlhzE20tzYsS1qKYaQ4K3iISQJmOJAPDk8yW7JZV0HNskI3trWzf2NaAM5MgSyfjC9Z5a9pcRKRuMp1JsjnHmZHJeceL97uv1v1uWSyViM0feU8pYU1EpG66O/NZ5Aunzk8OTTBwYVJT5lJWOhFnZEL3vEVEGqKw5enCjPNnCvuZK1lNykgn44xO5Ufezjllm4uI1FMxeC8Yef/y+BCd6+Ls2NTeiNOSgEsnYsWR9/j0LDkXza1RQcFbRAIo2dLMhraWRdXFeo4PctPV62lq0v1uWSydnKvpXdzXXCPvxcxsvZk9bma93vfOJZ6bNrM+M/taNa8pImtDvjTo3Mi7f3iCk0MTmjKXilKJGCMTMzjnihXFUhGsKAbVj7zvB550zu0AnvT+XMlfA/9d5euJyBrRvWCt91z9biWrSXnpRJycg4vTs8XNWlKaNi/rTmC/93g/sKfck8zsPUAX8FiVrycia0SmI0nfWxM4l1/r3XNsiFQixts3pxt8ZhJUpcVJNG2+tC7n3ACA933TwieYWRPwFWBvla8lImtIpiPJVDbH+bFpIL+z2k1XradZ97ulgsL+5qOT2bmRd0SD97JXZWZPAJvL/OgLK3yNe4FHnXMnl9tUwczuAe4B2LZt2wqbF5EoypSs9XbOcfz8RX7/pisafFYSZOlkobLYzFwt74hOmy97Vc65j1f6mZmdMbMtzrkBM9sCnC3ztPcBHzSze4F2oMXMxpxzi+6PO+ceBB4E2L179+J9EUVkzejunFvrfWIon3WuZDVZSmllsZGIJ6xV+5HkR8BdwAPe94cXPsE595nCYzO7G9hdLnCLiJTKFIL38DgnhsZpb43xjq263y2VzdX0LrnnHdGRd7X3vB8AbjWzXuBW78+Y2W4z+2a1Jycia1c6ESeViNH31gQ9x4Z4z5WdxJq1NYVUNpewlmVsMsu6lubI5khU9ZHEOTcIfKzM8eeAz5U5/i3gW9W8poisHZmOJIdPXaD37BifujHT6NORgJur6T0T6YpioB3WRCTAujuTPH9yGND9bllea6yZRLyJkcmsV1FMwVtEpO6msrni4z/53q84cKivgWcjYZCvLDbD6FSW9ogmq4GCt4gE1IFDfTz9+mDxz/0XJvn8Qy8qgMuSCjW9RydnSGvaXESkvvYdPEo2N3/F6MTMLPsOHm3QGUkYpJPxYsKaps1FROqsf0E50OWOi0B+2ryQsKbgLSJSZ1u9mt4rPS4ChbKg+YS1qG7QAgreIhJQe2/bSTLePO9YMt7M3tt2NuiMJAzSiRjD49P5bPMI3/OO7pWJSKjt2ZVf173v4FH6hyfY2pFk7207i8dFykkl4rw1nt8aNcoJa9G9MhEJvT27MgrWsiqF4iQQ3a1RQdPmIiISIemS+9xRnjZX8BYRkcgo7G8O0a0oBgreIiISIaX3uTVtLiIiEgLzR94K3iIiIoFXOvJW8BYREQmBeQlrmjYXEREJvsK0uRm0tUQ3eEf3ykREZM1pjTXR0txEa6yJpiZr9OnUjEbeIiISGQ8/3082l2N0Ksv7H/hpZEvIKniLiEgkHDjUx+cfepFCJdm+4YnI1oBX8BYRkUjYd/AoEzOz845FtQa8greIiETCWqoBr+AtIiKRsJZqwCt4i4hIJKylGvBaKiYiIpGwlmrAK3iLiEhkrJUa8Jo2FxERCRkFbxERkZBR8BYREQkZBW8REZGQUfAWEREJGQVvERGRkFHwFhERCRkFbxERkZBRC+ChVAAABL1JREFU8BYREQkZBW8REZGQMedco8+hLDM7B7zZ6PO4RBuB840+iYBS31SmvqlMfVOe+qWysPbNlc65y5d7UmCDd5iZ2XPOud2NPo8gUt9Upr6pTH1Tnvqlsqj3jabNRUREQkbBW0REJGQUvGvjwUafQICpbypT31SmvilP/VJZpPtG97xFRERCRiNvERGRkFHwXiEzu8LMfmZmR8zsJTP7U+/4ejN73Mx6ve+d3nEzs38ys9fM7AUzu9E7foOZPe218YKZfbqR11Utv/qlpL20mfWZ2dcacT1+8rNvzGybmT3mtfWymV3VmKvyh8998yWvjSPec6xR1+WHS+ibt3vvKVNmdt+Ctj5hZke9fru/EdfjF7/6pVI7oeOc09cKvoAtwI3e4xTwKnAd8CXgfu/4/cA/eI9vB34CGHAL0OMdvxbY4T3eCgwAHY2+vkb3S0l7/wh8D/hao68tSH0D/Bdwq/e4HVjX6OsLQt8Avw78Amj2vp4GPtLo66tz32wC3gv8LXBfSTvNwOvAdqAFOAxc1+jrC0C/lG2n0de32i+NvFfIOTfgnPuV93gUOAJkgDuB/d7T9gN7vMd3At92ec8AHWa2xTn3qnOu12unHzgLLLsgP6j86hcAM3sP0AU8VsdLqBm/+sbMrgNizrnHvbbGnHPj9bwWv/n4e+OABPng1ArEgTN1u5AaWG3fOOfOOueeBWYWNHUT8Jpz7phzbhr4d6+NUPKrX5ZoJ1QUvC+BN2W5C+gBupxzA5D/pSD/aQ/yvwwnS/7aKRb8gpjZTeTfdF6v7RnXRzX9YmZNwFeAvfU633qq8nfmWmDYzB4ys0Nmts/Mmut17rVWTd84554GfkZ+BmsAOOicO1KfM6+9FfZNJcu+B4VVlf1SqZ1QUfBeJTNrB/4D+DPn3MhSTy1zrJja740avgP8oXMu5+9Z1p8P/XIv8Khz7mSZn4eaD30TAz4I3Ed+GnA7cLfPp9kQ1faNmV0D/BrQTT4wfdTMPuT/mdbfKvqmYhNljoV+eZEP/eJrO42i4L0KZhYn/4/9b865h7zDZ0qmfbeQnwaH/KfcK0r+ejfQ7z0vDfwY+KI3BRhqPvXL+4A/NrM3gC8DnzWzB+pw+jXlU9+cAg55059Z4AAwL9EvjHzqm08Bz3i3EsbI3xe/pR7nX0ur7JtKKr4HhZVP/VKpnVBR8F4hL4P1X4AjzrmvlvzoR8Bd3uO7gIdLjn/Wy5K9BbjgnBswsxbgh+Tv3/2gTqdfM371i3PuM865bc65q8iPML/tnAt7dqwvfQM8C3SaWSE34qPAyzW/gBrysW9OAB82s5j3hvxh8vcwQ+sS+qaSZ4EdZna1977ze14boeRXvyzRTrjUOiMuKl/AB8hPOb0APO993Q5sAJ4Eer3v673nG/B18vezXwR2e8f/gHwCxfMlXzc0+voa3S8L2rybaGSb+9Y3wK1eOy8C3wJaGn19Qegb8hnV3yAfsF8Gvtroa2tA32wmP8oeAYa9x2nvZ7eTz6Z+HfhCo68tCP1SqZ1GX99qv7TDmoiISMho2lxERCRkFLxFRERCRsFbREQkZBS8RUREQkbBW0REJGQUvEVEREJGwVtERCRkFLxFRERC5v8B4p92oYj+RQUAAAAASUVORK5CYII=\n",
      "text/plain": [
       "<Figure size 576x432 with 1 Axes>"
      ]
     },
     "metadata": {
      "needs_background": "light"
     },
     "output_type": "display_data"
    }
   ],
   "source": [
    "wine = alcohol_consumption.TotalWine\n",
    "wine_diff = wine.diff(4).dropna()\n",
    "\n",
    "plt.figure(figsize=(8, 6))\n",
    "plt.plot(wine_diff, \"-o\")\n",
    "plt.title(\"Does this look stationary ?\");"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "It looks stationary enough. Let's call the result `time_series`."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 3,
   "metadata": {},
   "outputs": [],
   "source": [
    "time_series = wine_diff"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## Autoregression model\n",
    "\n",
    "In an autoregression model, values are modelled as a linear combination of the $p$ past values. An autoregressive model of order $p$, that is usually indicated as $AR(p)$ model, can be written as:\n",
    "\n",
    "$$\n",
    "y_{t} = c + (\\phi_{1}y_{t-1} + \\phi_{2}y_{t-2} + \\dots + \\phi_{p}y_{t-p}) + e_{t},\n",
    "$$\n",
    "\n",
    "where\n",
    "\n",
    "* `c` is the mean of the time-series\n",
    "* `e_t` is the noise\n",
    "\n",
    "The library `statsmodels` contains a `tsa` module which implements the autoregressive model.\n",
    "\n",
    "**Task**:\n",
    "\n",
    "* Define an `AR` model using `tsa.AR` \n",
    "* select the order using the `select_order` method. You will need to specify a maximum `p` to consider and a criterion for deciding which model is \"best\" \n",
    "  * AIC and BIC are common parameters which encourage \"good fit\" while penalising having too many parameters (complex model). "
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 6,
   "metadata": {},
   "outputs": [
    {
     "name": "stderr",
     "output_type": "stream",
     "text": [
      "/anaconda3/envs/cspark/lib/python3.7/site-packages/statsmodels/tsa/base/tsa_model.py:225: ValueWarning: A date index has been provided, but it has no associated frequency information and so will be ignored when e.g. forecasting.\n",
      "  ' ignored when e.g. forecasting.', ValueWarning)\n"
     ]
    },
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "Optimal p: 5\n"
     ]
    }
   ],
   "source": [
    "#com+ add your code here\n",
    "ar = tsa.AR(time_series)\n",
    "optlag = ar.select_order(10, ic=\"aic\")\n",
    "print(\"Optimal p: \" + str(optlag))"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "### Having a look at the AR model\n",
    "\n",
    "Let's see if the AR(p) does a good job compared to the original time series. \n",
    "\n",
    "* use the `fit` method specifying the optimal lag found above\n",
    "* use the `predict` method to generate values starting at the optimal lag\n",
    "* plot the predicted results and the data, how does it look?"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 7,
   "metadata": {
    "scrolled": true
   },
   "outputs": [
    {
     "data": {
      "image/png": "iVBORw0KGgoAAAANSUhEUgAAAs4AAAD8CAYAAABnwSKjAAAABHNCSVQICAgIfAhkiAAAAAlwSFlzAAALEgAACxIB0t1+/AAAADl0RVh0U29mdHdhcmUAbWF0cGxvdGxpYiB2ZXJzaW9uIDMuMC4xLCBodHRwOi8vbWF0cGxvdGxpYi5vcmcvDW2N/gAAIABJREFUeJzs3Xl4VOX1wPHvO5PJvu9kI6yRnUBQFlnEBRFBRMSldWlFrda60x9Wa9G2arWtS611QetaFQVRREVBQUV2wo6sISEJCUnIvs5yf3/cTMgySSbJTBLC+TyPT8ydO/e+JJmZc9973nOUpmkIIYQQQgghWmbo6gEIIYQQQghxJpDAWQghhBBCCCdI4CyEEEIIIYQTJHAWQgghhBDCCRI4CyGEEEII4QQJnIUQQgghhHCCBM5CCCGEEEI4QQJnIYQQQgghnCCBsxBCCCGEEE7w6OoBNCc8PFxLTEzs6mEIIYQQQogebtu2bfmapkW0tl+3DZwTExPZunVrVw9DCCGEEEL0cEqpdGf2k1QNIYQQQgghnCCBsxBCCCGEEE6QwFkIIYQQQggnSOAshBBCCCGEE7rt4kAhhBBCCNFQSUkJJ0+exGw2d/VQzjh+fn7ExcVhMLR/3lgCZyGEEKITLU/N4plVB8guqiQm2IcF05KYnRzb1cMSZ4CSkhJyc3OJjY3Fx8cHpVRXD+mMYbPZyMrKIj8/n8jIyHYfxyWpGkqpS5VSB5RSh5VSCx08nqCU+k4plaqU2qWUuswV5xVCCCHOJMtTs3ho2W6yiirRgKyiSh5atpvlqVldPTRxBjh58iSxsbH4+vpK0NxGBoOBqKgoiouLO3acjg5EKWUE/g1MBwYD1ymlBjfa7RFgiaZpycC1wEsdPa8QQghxpnlm1QEqzdYG2yrNVp5ZdaCLRiTOJGazGR8fn64exhnLZDJhsVg6dAxXzDifCxzWNO2opmk1wAfAFY320YDA2v8PArJdcF4hhBDijJJdVNmm7UI0JjPN7eeKn50rAudY4Hi97zNrt9W3CPilUioT+AL4nQvOK4QQQpxRYoIdzxY2t10I0b24InB2FL5rjb6/DnhT07Q44DLgHaVUk3MrpW5TSm1VSm3Ny8tzwdCEEEKI7mPBtCR8TMYG23xMRhZMS+qiEQnhOomJiaxevbqrh+FWrgicM4H4et/H0TQV4xZgCYCmaRsAbyC88YE0TXtV07QUTdNSIiIiXDA0IYQQovuYnRzLk3OGYay9ZRzkY+LJOcOkqobo8TqaW9xduCJw3gIMUEr1UUp5oi/++6zRPhnAhQBKqUHogbNMKQshhDjrzBoRg72M7DVj4iVoFp1ieWoWE576lj4LVzLhqW9dXsnlhhtuICMjg5kzZ+Lv78/TTz+NUorXX3+dhIQEpk6dytq1a4mLi2vwvPqz1Dabjaeeeop+/foRFhbGvHnzOHXqlEvH2VEdDpw1TbMAdwGrgP3o1TP2KqUeV0rNqt3tAeBWpdRO4H3gZk3TGqdzCCGEED1eYUUNZqv+EZhRUNHFoxFng84og/jOO++QkJDAihUrKCsrY968eQCsW7eO/fv3s2rVqlaP8cILL7B8+XLWrVtHdnY2ISEh/Pa3v3XZGF3BJQ1QNE37An3RX/1tj9b7/33ABFecSwghhDiT5ZRUAWA0KDJOSeAs2u+xFXvZl13S6n6pGUXUWG0NtlWarfz+4128vzmjxecOjgnkTzOHtHuMixYtws/Pz6l9X3nlFV588cW6WelFixaRkJDAO++8g4dH9+jZ1z1GIYQQQpwlcmsD5yExgaTllaNpmpQYE27VOGhubbsrxcfHt75TrfT0dK688soGLbGNRmNdt8TuQAJnIYQQohPlFFcDkNI7lF2ZxRRVmAnx8+ziUYkzkbMzwROe+pYsB7XCY4N9+PD2cS4bj6MLwPrb/Pz8qKg4fZfFarVSv4pafHw8b7zxBhMmdN8kBZe03BZCCCGEc3JLqlAKUhJDACRdQ7hdZ5VBjIqK4ujRo80+PnDgQKqqqli5ciVms5m//OUvVFdX1z3+m9/8hocffpj09HQA8vLy+PTTT106xo6SwFkIIYToRLklVYT5edE3Qs/7lMBZuJu9DGJssA8KfabZHWUQH3roIf7yl78QHBzMxx9/3OTxoKAgXnrpJebPn09sbCx+fn4Nqmzcc889zJo1i0suuYSAgADGjh3Lpk2bXDrGjlLdtbhFSkqKtnXr1q4ehhBCCOFSN/93M/ll1Xx42ziG/GkVC6Yl8dsL+nf1sMQZYP/+/QwaNKirh3FGa+5nqJTapmlaSmvPlxlnIYQQohPlFFcRFeCNn5cH4f6eHJcZZyHOGBI4CyGEEJ3oZGk1UUHeAMSH+kqqhhBnEAmchRBCiE5SbbFyqryG6EA9cE6QwFmIM4oEzkIIIUQnOVmiVxCoHzhnF1Vi7oR6ukKIjpPAWQghhOgk9uYnkYFegJ6qYdMg20GNXSFE9yOBsxBCCNFJ7O22o4NOzziDlKQT4kwhgbMQQgjRSXKKawPnQAmchTgTSeAshBBCdJLckio8PQwE+ZgAiAr0xtNokMBZiDOEBM5CCCFEJ8ktqSY60BulFABGgyIuxEdqOQvRgptvvplHHnnEqX0TExNZvXq128YigbMQQgjRSXJKqurSNOyklrMQZw4JnIUQQohOkltSVdf8xC4h1JeMAgmchTgTSOAshBBCdAJN02rbbXs12J4Q6ktJlYXiCnMXjUycFXYtgWeHwqJg/euuJS4/RWJiIs888wzDhw/Hz8+PW265hdzcXKZPn05AQAAXXXQRhYWFAHz22WcMGTKE4OBgpkyZwv79++uOk5qayqhRowgICOCaa66hqqqqwXk+//xzRo4cSXBwMOPHj2fXrl0u/7c0RwJnIYQQohOUVFqottjqStHZxUtlDeFuu5bAiruh+Dig6V9X3O2W4Hnp0qV88803HDx4kBUrVjB9+nSeeOIJ8vPzsdlsvPDCCxw8eJDrrruO5557jry8PC677DJmzpxJTU0NNTU1zJ49mxtuuIFTp05x9dVXs3Tp0rrjb9++nV//+te88sorFBQUcPvttzNr1iyqq6td/m9xxKNTziKEEEKc5ew1nKMCm6ZqgB44D4sL6vRxiTPYlwshZ3fr+2VuAWujwNJcCZ/eBdveavm50cNg+lNOD+l3v/sdUVFRAEycOJHIyEiSk5MBuPLKK1mzZg0Gg4EZM2Zw8cUXA/Dggw/y/PPP89NPP2EwGDCbzdx7770opZg7dy7//Oc/647/2muvcfvtt3PeeecBcNNNN/HEE0+wceNGJk+e7PQ420sCZyGEEKITNBc4x4f6AHC8UGachZs0Dppb294B9qAZwMfHp8n3ZWVlZGdn07t377rtBoOB+Ph4srKyMBqNxMbG1lWeARrsm56ezltvvcW//vWvum01NTVkZ2e7/N/iiATOQgghRCewt9tuXFUjwNtEqJ+npGqItnN2JvjZobVpGo0ExcOvVrp2TE6IiYlh9+7TM+WapnH8+PG6gDkrKwtN0+qC54yMDPr16wdAfHw8Dz/8MA8//HCnjxtclOOslLpUKXVAKXVYKbWwmX3mKaX2KaX2KqX+54rzCiGEEGeK3NqugZGBXk0eiw/1lVrOwn0ufBRMPg23mXz07V1g3rx5rFy5kjVr1mA2m/nHP/6Bl5cX48ePZ9y4cXh4ePDCCy9gsVhYtmwZmzdvrnvurbfeyssvv8ymTZvQNI3y8nJWrlxJaWlpp4y9w4GzUsoI/BuYDgwGrlNKDW60zwDgIWCCpmlDgHs7el4hhBDiTJJTUkWIrwlvk7HJYwlSy1m40/B5MPMFfYYZpX+d+YK+vQskJSXx7rvv8rvf/Y7w8HBWrFjBihUr8PT0xNPTk2XLlvHmm28SEhLChx9+yJw5c+qem5KSwmuvvcZdd91FSEgI/fv358033+y0sStN0zp2AKXGAYs0TZtW+/1DAJqmPVlvn6eBg5qmLXb2uCkpKdrWrVs7NDYhhBCiu5j/1hYyCyv56t5JTR57+qufefX7o/z850vxMErBK+HY/v37GTRoUFcP44zW3M9QKbVN07SU1p7vildnLFA/cSazdlt9A4GBSqn1SqmNSqlLHR1IKXWbUmqrUmprXl6eC4YmhBBCdA+5JdVNFgbaJYT6YrFpnCiucvi4EKJ7cEXgrBxsazyN7QEMAKYA1wGLlVLBTZ6kaa9qmpaiaVpKRESEC4YmhBBCdA+O2m3b2UvSSZ6zEN2bKwLnTCC+3vdxQOOaIJnAp5qmmTVNSwMOoAfSQgghRI9nttrIL6tu0m7bTpqgCHFmcEXgvAUYoJTqo5TyBK4FPmu0z3LgAgClVDh66sZRF5xbCCGE6PbySqvRNIhyUFEDoFeQNx4GJYGzEN1chwNnTdMswF3AKmA/sETTtL1KqceVUrNqd1sFFCil9gHfAQs0TSvo6LmFEEKIM0FzNZztPIwGYkN8JHAWoptzSQMUTdO+AL5otO3Rev+vAffX/ieEEEKcVXKb6RpYX4LUchZOsNlsGAxSeaU9OlpJDlzUAEUIIYQQzcspbj1wjpdazqIVfn5+ZGVlUVNT45Ig8GyiaRoFBQV4ezf/GnSGtNwWQggh3Cy3tBqTURHm59nsPgmhvhRWmCmpMhPoberE0YkzRVxcHPn5+aSnp2OxWLp6OGccb29v4uLiOnQMCZyFEEIIN8striIywBuDwVEFV139knRDYoI6a2jiDGIwGIiMjCQyMrKrh3LWklQNIYQQws1ySqqarahhJ7Wchej+JHAWQggh3EwPnFvOrZRazkJ0fxI4CyGEEG52soV223ZBPiaCfEwSOAvRjUngLIQQQrhRWbWFsmoL0c10DawvIdSXjFOVnTAqIUR7SOAshBBCuJG9FF1zzU/qk1rOQnRvEjgLIYQQbmRvfhLZyuJA0POcMwsrsNqkRq8Q3ZEEzkIIIYQbtdZuu76EUF/MVo2c2ucIIboXCZyFEEIIN8pxot22nb0kXUaBpGsI0R1J4CyEEEK4UW5xFQFeHvh5td5zTGo5C9G9SeAshBBCuFFuSTVRTlTUAOgV7I3RoKQknRDdlATOQgghhBvllFQ5ld8MYDIaiAn25nihBM5CdEcSOAshhBBulOtE18D69FrOEjgL0R1J4CyEEEK4idWmcbK0mignStHZSS1nIbovCZyFEEIINykor8Zq05zqGmgXH+pLflkN5dUWN45MCNEeEjgLIYQQbpJbXA04V4rOLj6ktrKG5DkL0e1I4CyEEEK4SU4bmp/YSS1nIbovCZyFEEIIN2lL8xO7usBZ8pyF6HYkcBZCCCHc5GRJFQYF4f6eTj8n2NdEgJeHLBAUohtySeCslLpUKXVAKXVYKbWwhf3mKqU0pVSKK84rhBBCdGc5xVVEBHjhYXT+41YpRbyUpBOiW+pw4KyUMgL/BqYDg4HrlFKDHewXANwNbOroOYUQQogzQU4bazjbSS1nIbonV8w4nwsc1jTtqKZpNcAHwBUO9vsz8DRQ5YJzCiGEEN3eyZLq9gXOYb4cL6zEZtPcMCrhjOWpWUx46lv6LFzJhKe+ZXlqVlcPSXQDrgicY4Hj9b7PrN1WRymVDMRrmvZ5SwdSSt2mlNqqlNqal5fngqEJIYQQXact7bbriw/1pcZi42RptRtGJVqzPDWLh5btJquoEg3IKqrkoWW7JXgWLgmclYNtdZfISikD8CzwQGsH0jTtVU3TUjRNS4mIiHDB0IQQQoiuUWW2UlxpblPzEzuprNG1nll1gEqztcG2SrOVZ1Yd6KIRie7CwwXHyATi630fB2TX+z4AGAqsVUoBRAOfKaVmaZq21QXnFy1YnprFM6sOkF1USUywDwumJTE7Obb1JwohhOiQnGI9MzEywPl223b1A+dz+4S6dFyiddlFlW3aLs4erphx3gIMUEr1UUp5AtcCn9kf1DStWNO0cE3TEjVNSwQ2AhI0dwK51SSEEF0n1978pB0zzrHBPiglM85dJSbYp03bxdmjw4GzpmkW4C5gFbAfWKJp2l6l1ONKqVkdPb5oP7nVJIQQXafZroG7lsCzQ2FRsP5115Imz/X0MBAT5CO1nLvIgmlJeBobZqL6mIwsmJbURSMS3YUrUjXQNO0L4ItG2x5tZt8prjinaJ3cahJCiK5jn3GOrB8471oCK+4Gc+37cPFx/XuA4fMaPD8+1EdmnLvI7ORY3lh/lN2ZJWhAkI8Hj80aKqmOQjoH9mRyq0kIIbpOTnE1PiYjgd715qjWPH46aLYzV+rbG+nptZy7c7m3kiozP+eUcdP4RAK8PLhiZKwEzQKQwLlHWzAtCZPcahJCiC6RW1pFdJA3tQvjdcWZjnd2sD0h1Je80moqa6wOnnBm6+5rcL7anUONxcbs5FgSwnxJL+i5FzCibSRw7sFmJ8cyrm9Y3fc+JgNPzhkmV81CCNEJcouriApsVFEjKM7xzn7hTTbF11bWOF7Y84K27r4G55PULPqE+zEiLojeYT175l+0jQTOPVyl2cqohGCuHh2Hp4eRy4f36uohCSHEWcFhu+0LHwVlbLSngvI8+P7vYLPVba0rSdcDZzu78xqcnOIqNqYVcMXIGJRSJIT6kVlYgVW6OAokcO7RrDaNvdklDIsNYkpSJMWVZnZmFnX1sIQQosfTNI2TJdVNK2okng+aDTz9AQVB8TDrBRg6F779M7x/DVScAnp2E5TuvAbns51ZaBrMHqnfne0d5ovZqnWLoF50PQmce7C0/DIqaqwMjQ3i/P7hGBSsOyCtzIUQwt0KK8zUWG1NZ5w3vQJKwR3rYVER3LcHRt0IVy2GGf+Ao2vhlUmQuY1QP0/8PI09MnBeMC0JH1PDmffusgZneWo2I+KDSQz3A6B3D76AEW0ngXMPtjurGIBhcUEE+ZoYlRDC2oMSOAshhLvZuwY2aH5SXQbb/guDZkJIYsMnKAVj5sOvvwIUvDENtWUx8SE+ZPbAHOfZybE8ceVQ7MsmDQqeuLLry70dzC1l34kSZo+MqduWEKYHzrJAUIAEzj3a7swSvE0G+kf4AzB5YAS7MovJL6vu4pEJIUTPZq/h3GBxYOq7UFUM437X/BNjR8Pt66DfVPjiQR63PEteQb6bR9s1xvcPRwOSogKwaTC6d9e3Fl+emoXRoLh8eExdo5rY52NY73U3fgeXdvXwRDcggXMPtiermEG9AvEw6r/mKUmRAPxwSGadhRDCnU4HzrUzzjYrbHwJ4s+D+DEtP9k3FK77AC78Eynl63i2+H60H55ttdvgmeZoXjkA88bEA7D52KmuHA42m8anO7I5v384EWmf6o1pio+j0IhV+Vx69Mke8XNvMyc6XZ5NJHDuoWw2jb3ZxQyLDarbNiQmkHB/T9ZKnrMQQriVvd12ZEBt4Pzz51CUDuPucu4ABgNMvJ9vxrxKOIWwZpHeZRDtdLfBMzyAScvXA+dLBkcR5GNiS1rXBs7bMgrJKqpkdnKMw0Y1Xlq1w0Y1PZq902UP+9vrCAmce6ij+eWU1y4MtDMYFJMGRPD9wTwpqyOEEG6UW1JFmJ8nnh61H7M/vajnNZ8zo03H8ew3mXLNG9X4gWa6DZ5J0vLL8PIwEBvsw5jEELakd23g/ElqFj4mI5cMjm62UY3WXAObnqoNnS7PFhI491B77AsD6wXOAJOTIiisMNctHBRCCOF6uSXVp9M0jm+GzM0w9k4wNK7h3LL4UF+iVKHjB8/wIC4tv5zEMD8MBkVKYihH88q7bA1OjcXGF7tPcMmQKPy8PJptVGMLPMsaiLWh0+XZQgLnHmp3VjFeHgYGRPo32D5xQARKytIJIYRb5RRXna6o8dO/wDsIRv6izceJC/EhW2vaVRCAgKgOjLDrHc0vp09tybcxifrCwK1dlOe87mAeRRXmutrNjL2jyT4VmicZyQ928si6WHOdLpvbfhaQwLmH2p1ZzOCY0wsD7UL9PBkeF8zagye7aGRCCNHz5dq7Bp5K0/ObU34NXv6tP7ERb5ORV02/pEZ5NX2wsliv+3wGslhtZBRU0CdCD5yHxQbhbTKwOa2Z2XU3W74ji1A/T84fUHuRcnwzGEwQoJelsxk8WWiez87gS7pkfF3mwkfB1KgpjdFT336WksC5B3K0MLC+KQMj2Hm8iMLymk4emRBC9HzVFisF5TV6KbpNL+stts+9vd3H+zniUv4TcI/eZdDebXDak3rO9DtzYOsbLht7Z8ksrMRi0+pmnD09DIyMD2ZLF8w4l1aZWb0vl8uH98JkNEDGJti3HCbeDw/sh/N+gzJ68Llt/NlXy3n4PJj2RL0NCnol69vPUhI490COFgbWNyUpApsGPxzumbVBhRCiK+WV6nm68T7VsP0dGDYXAnu1+3jxob68XzVW7zJo7zY47k645Wu93vPn98FXD+kl784QaQV6RY2+tYEzwLmJoezNLqas2tKpY/lqTw7VFhtXjIwFTYOvHwb/aBh/t75DTDLKXMF5/vmknyrv1LF1C8EJ+tcbP9O7XObugZqz7AKiHgmce6DmFgbaDY8LJsTXJHnOQgjhBvYaziNPLgdzOYz7bYeOlxDqS05JFVXmRoGxd6Be7/m8O/Qa0e9fB9Wl7T9RJ9brTaut4dynXuCckhiKTYPUjM5N1/h0RzYJob6MSgiGvcsgcwtMfeR0ak1MMgDn+x0n42ybcQbITtW/9hoBw67W/6YPftm1Y+pCEjj3QM0tDLQzGhQTB0Sw7mAeNheUpVuemsWEp76lz8KVTHjqW5anZnX4mEKItpPXYveQU1yNCQsJh96BPpMheliHjpcQpueYZhZWNn3Q6AHTn4IZ/4TDq+H1aVCU0faTdHK93rT8cgK9PQj186zbNqp3CAZFp9ZzPllSxU9H8pk9MgZlqYbViyBqKIy8/vROYf3B058RxjTST52NgfMOCO0HPsHQe7ye9737464eVZeRwLkH2t2oY6AjkwdGkF9Wzb4TJR061/LULB5atpusoko0IKuokoeW7ZYP7O5Iuj/1aPJa7D5ySqq43LABU0UujG+hvbaTEkJ9ATjeUtA25hb45cd6mbDXpsLxLW07SSfX603LL6dPhD9Kna5Q7e/lwZCYoE7tIPjZzmxsGlyRHAubX9EvOi75S8OygQYj9BpBf/Mh8kqrqajp3FSSLpedWjfrjsEIQ+fAoW+gsmsWcnY1CZx7GJtNY192SbNpGnaTBkYAegmejnhm1QEqG90+rDRbeWbVAecOIMFc55DuTz1eh1+LwmVOFldyq+lLtIhzoP9FHT5efG3gnNHabGe/qTB/NXj6w5sz9FnBlt5jy/Jg10fwyR217w0OuKleb1p+OX3CfJtsH5MYSmpGETUWm1vO29jyHVkMiw2in281fP8PGHAJ9Lug6Y4xyUSUH8SIleOnHMz891SluVCSdTpwBj1n32aGfZ913bi6kATOPUxaQTll1ZZWA+eIAC+Gxgay9kDHytJlFzl+A2luewMSzHUe6f7U43XotShcKiBnA4PVMdTYO0E16fnXZhH+XnibDK0HzgARA2H+GogdDUtvgeV3NHyP/fQueO8aeHki/L0/LJuv56s2Ljlm54Z6vVVmK1lFlfQJb5pOeG6fEKottk5p0nX4ZBl7skqYnRwL656CmjK4+M+Od45JxmirZoDKIr3gLFogeGKH/rV+4NxrpJ66sefsTNdwSeCslLpUKXVAKXVYKbXQweP3K6X2KaV2KaXWKKV6u+K8oin7wsDmKmrUN2VgJNsziiiuNLf7fDHBjt9sjQZVN5ZmSTDXeaT7U49WbbGebu3cSHiAg/q/HSS51C0bf/J9igzBMPwalxxPKUVCqK9zgTOAXxjcuBxMvmBrlFZgrYZDX4FXIEz9I9z6HSw4AjNfaBo8m3zcUq/3WG3gaa/hXF9KbSOUzihL9+mOLAwKZseVwZbXYfRNEHmO451rA8dhhqPO/x7cpFNff9mp6CXohp/eppS+SDDtByg54b5zd1MdDpyVUkbg38B0YDBwnVJqcKPdUoEUTdOGAx8DT3f0vMKx3ZnFeHoYGBDVeqH9yUkRWG0a6ztQlm7+xD5NtnkaDfh6Gpj97/W8sOYQFmszt9wkmOs80v2px7JYbdzz/g6qLTZMxoazmwooqaxxaTc2yaVuRd4BRlVv4ceQ2WDydtlhE0J9W85xbszDC63xxEQtDQW/WgmTHoTYUXre6vB5evAcFF/7fB/9ezfU6z2W37QUnV24vxd9w/3c3kFQ0zSW78hiQv9wwn76q36RMeUPzT8hpA94BZFiOtaltZw7/fWXnQrhA8EroOH2YXMBDfZ+4p7zdmOumHE+FzisadpRTdNqgA+AK+rvoGnad5qm2f/SNgJn7ae1u68U7QsDTS0sDLRLjg8mwNujQ2XpMgsrUUBUoBcKiA324em5w/n+91O5bFgv/vnNQa56eQNH8sqaPlmCuc5z4aNg8Gi0UcGUh7pkOMI1bDaNBR/v4qu9Ofzx8sE8M3cEscE+da/FR2cOJibYlxte38z3HVzPYCe51C3TfnqRKs3EgTjXBpzxtTPOmuZ8JaRcHLfqbm47w+fpNaIHXKJXknBTk4ujtYFzooPAGfQ85y3HCl1S9ak52zMKOX6qkvlxmXqqysT7wT+i+ScYDBAzgmTjsS6trNGprz9Ng+xUMryTmsYt4QP08nS7P3L9ebs5VwTOsUD9VQWZtduacwvgsACgUuo2pdRWpdTWvLyeV2PY3VeKesfAEoY7kaYB4GE0MHFAOOsO5rXpzdiuuNLMB5szuGJkDJv+cBFpT81g/cKpzE6OJdjXkxeuS+bF65NJLyjnsud/4L/r0xq+ETqqbeqmW4NnvYGX6t3LTL6AAt9wQIPs7V09MtFOmqbxx0/38ElqFg9eMpBbzu/D7ORY1i+cWvda/NWEPiy5fRyJ4X7c8tYWvtzd8duqkkvdgrI82PUhS62TCAxrf8MTRxJCfamo0TsSOkPTNJ6ouZoKzbPB9grNkydrrm75yaF94dRRPXByg7S8ciIDvPD3anwxrxvTJ5TiSjOHTjqYcHGR5anZ+Jjg/KPP6rPsY+9s/Um9RtLHlkZ2fpHbxtWaTn39lZ6AslzeyQh1HLcMu1r/DCk44vpzd2PqJswMAAAgAElEQVSuCJwdrXxw+GpTSv0SSAGecfS4pmmvapqWomlaSkREC1d+Zyh3Xyk6uzCwvikDI8kpqeJAbtuL5r+3KZ3yGiu3TerX7D6XD4/h63snMb5fGI+t2McvX99Elv0FXl6bIuIfrX/1CnDbrcGz3s73wVrNLWoRfareY4LtNQ73uwm2LIY9y7p6dKKNNE3jiS/2896mDO6Y0o/fXtC/2X0jArz44NaxDIsN4rf/285HW5upnuDEOT/dkYVSMMvwIz963s1Rr+v50fNuZhl+bHa9w1lly2KUtZo3rJcSGeja3PIEZytroK91+cXiTXxmO5+F5vlk2sKxaYpMWzgLzfPZGnhxywcI7ac3uSjLdcXQm0jLL2/Q+KSxc2vznN1Vls5stbFy9wn+ELMTY+5uuPBPzqXVxCRj0sz4FbeQguhmzb3O3PL6q218ss2c2GBzXdwyZA6gYM9S15+7G3NF4JwJxNf7Pg7IbryTUuoi4GFglqZp1S447xnH3VeKbVkYaDc5Sb9AWdvGdI1qi5X/rj/GxAHhDI4JbHHfyEBv3rh5DE/OGcbO40Vc+uz3LN1yFG372zBwOjx4ABLG6XlUEjS7ns1G2ff/Zoc2gDUlcXWzBlcevIRTISPgs7u7zYyBLDpzzvNrDvHaD2ncNK43v5+W1KAWriNBvibenX8eE/qHs+DjXbzxY1qbzrcrs4i5L2/gng92cKPfZp4yLSbOkI9BQZwhn7+ZFvPc4EMd+Sed+cyVsOU1CmIv4IgWS3Sg6/KbwblaztlFldy/ZAczX/yR/SdKmJMcw9fGSZxf8wJ9q9/j/JoX+MY4mQXTklo+WWhf/eupo64afgNp+eX0dbAw0C4+1IeoQC+3NUL54VAeFeWlXF3yJsSMgqFXOffE2gWCgznCieIq1w7KydKsC6Yl4WFo+Hr3MRlb/522R3YqFs3APq1pPYfsokoIioXeE/R0DTfdneiOXBE4bwEGKKX6KKU8gWuBBsX9lFLJwCvoQXPH6p+dwdx9pdiWhYF2UYHenBMd0OY8509Ts8krrea2SX2d2l8pxXXnJvDlPZMY1CuQ7z55A1V+kjUBM5nw1Le8eDQKS1Yqn28+2KZxCCcc+Rb/8nTeMF/SYHOpWTG/4rd657GPbgKziz8I2qg9qUydHWh3h8D+1e+P8NzqQ8wdHcefZg5pNWi28/X0YPFNKVw6JJrHP9/H86sPtZqidbK0igUf7WTWi+tJLyjn6auG8ye/pfiqhukCPqqGMUf+1e5/U4+w8wOoKGB3/A0ARAe5NnCOC6mdcXawMK20yszTX/3MBX9fy+e7TnDbpL6sXXAB/7wmmafmDKdX7Vh8PY08OWeYXn6tJWHuC5yLK8wUlNeQGNZ84KyUIiUxlC3HTrUrjbA1n6Rm8zufr/CuzIVpT+j5y84IScTsGcQwddS1CwTbUJr14sFReBoVXrVVdHxMTv5O2yM7lTRDAlU0vXtiNCi9sMCwqyD/IOTsdv35u6kOB86aplmAu4BVwH5giaZpe5VSjyulZtXu9gzgD3yklNqhlDorq2YvmJZEo0Xv+JgMLrtSbMvCwPqmJEWyNf0UZdXOdUOy2TRe/eEog3sFcn7/ZhaZNCMhzJf3bxvLHyLXk65FMn99IFlFlWy0DcIDG8tXLDsrZxndGZBpm17mpBbMl7bzmjyWWuwPV76iv+l91aSSZKd6etXPDlOZ/m/pLh7+ZDf/WXuEFTuzSc0oJK+0mk+2Z3bq6vLuUE3inY3pPPHFz8wY3ou/XTUcg6FtNYK9PIy8eH0yc0fH8ezqg/xl5X6HgUm1xcp/1h7hgmfWsnxHFrdP6st3D05h3ph4lFTDacpmgw3/hl4j2GPS22tHuXjG2cfTSGSAV4NUDbPVxjsbjjHlmbW8tPYIlw6NZs39k3lo+iCCfEwAzE6OZcNDF/KL8xKw2LS65lctCkrQFxO74U5Umr0UXQupGqCna5worjqd2uciZdUWduzbz3z1GQyaCb3HOf9kpbBGj2S4IY30Uy6s5dyG0qzvbEynwmzjg9vGMnFAOAOi/N0TNNcuDPSMH4WPydjgIU+jgUAfD36xeBMP7ElEM3icVYsEHWfmt5GmaV8AXzTa9mi9/+9466Qe4NKh0Sz4WOFlNFBRowcIt07q65I/evvCwNnJMW1+7uSBEby87gg/Hc7nkiHR+lXumsf1D8KgOH2xXr0Uim9/Psnhk2U8d81Ip2e76jPm7SemaDsvGG9Aq71222YbgFkzkqzt5ZlVB9zzRtBN2QMye9BoD8iADv8cCtL3EXb4G96zXIXZwcs9JtgHBk6FCffA+uch8fzaMkOdp7LGypKtx8kucjzjXW2x8cXuExRWtF5v3J57546/n7995Tiw76y/16XbMvnj8j1ceE4kz84bibGNQbOdh9HA01cNx9/Lg9d/TGNPVhHHCys5UVRFTLA304ZEs+bnk6QXVHDRoCgenjGoYZATFOe4y9zZXA3n8DdQcAjmLCbnaDVBPia8GwUbHbU8NYuiCjMfbcvkpyP5XDq0F2sPnORIXjnn9gnlvzMGMTwuuNnn3zgukfc2ZfDhluPcMaX5dSmAfhcqOMEtM85p+fqCv5ZSNUCvrAF6PWf7bLsrfL03hzu1JZiwwEWPtfn5nvGjSUr/kS/yCgEXtaRw8mK0osbCa98fZdLACJITQkgM8+OznU0yY100puNQUUDvYRN4MnkYC5fuospiIzbYhwXTkrh0aDSLfzjKi98dZoYaTsq2D/Gc8ijenib3jKcbkc6BnWjdwTzMVo1XbhjN3sem4edpJMdFeVLH2rEw0G507xD8vTxYezDPqVtGr35/lNhgH2YMb+eq8a2vg9GLNysm1G2qxJtdWl/OM/x81q3Od8eiUU3TWLY9k6/++zg1mpFTg67Hx9T05X7debXLE6b+EeLPgxX3QH7n5KoWltfw3OqDjH9qDX/6bC+ezdwpiQ32IfXRS9jz2DRW3TuJ129KYdHMxqXiT8sqqnTpLFVOcRV/++rnZnMaO+Pv9cvdJ1jw8U7G9wvj378Y1WyzE2cZDIo/zRzMtMFRbEorJLuoqnYWvYo31h+jqsbKO7ecy+KbUprODE5ZSOM14RaMaGdzNZyf/gWBsTBkNrkl1S7Pb7ZfXNfULkjLKqri9R/TKKky89qNKXx429gWg2aApOgAxvUN492N6VidKfMW2g9OuWHGOb8CgzrdRrw5SdEBBHh7sDmt0KXn37p5PfM81qHGzIewVi4gHDDEJmNSVmwn9rhuUE6WZn13YzoF5TXcc6G+GLh3mC/FlWaKKpyrtNImtQsDiUlmdnIsfSP8mXpOZF3lLG+TkbumDuDbB6aQETOdwOpc7v/7K3y154Rb0mu6EwmcO9GXu08Q7GtibN8w/Lw8uHx4DJ/vOkG5kykSLdndjoWBdp4eBsb3C2PdgTy0Vm4ZpWYUsvnYKX59fp82p4QAUF2q5wIOnYNPcFSDhzbZBjFCHaFvUMdb1Dbh5MKLrtBc4JVVVMmqvTnUWNq2eju3pIr5b23lj0s2cqVaS3XSFfz5lxfx5JzhdTV+o4O8CfQ2smxblv73ZzTB3P+C0RM+urnp34ALZRZWsOizvYx/6lueW32I0b1D+Pg343h67vAmtwTrL3rx9/IgKTqACwdFcfOEPsS2sDZgwlPfcu2rG1iy5TilVe3rjLknq5j7P9zB+X/7llfWHcHbwYUHuGk1ez3f/XySuz9IJTkhhNduTHHZTKZSij3ZJQ4fMxoVEwc0c0u/MB3QwC8CUJiNvijNSpbpLGsIW/895dgPkDAWjCZyS6qIcnF+s6OLawCT0cDFg6OcvvN30/jeZBVVsnq/E9UyQvvCqTSXL/pKyy8nLsQXL4+W/46NBkVK7xCXdhA8WVrFpdn/osboj5r8+/YdpHaBoP8pFwbOFzzcdFuj0qyVNVZe/f4o5/cPZ3RvfTa+d22euFsasmSngsEEUUMB/X07PqTpe11MsA83/+pOrEYfLuNHfvPudm54fTOHcku7xZoQd3BJqoZoXbXFypr9J7l0aHRdwDlvTBwfbj3Oyt0nmJcS38oRWrYnS18YODAqoPWdHZiSFMnX+3LBu+VbRq9+f5RAbw+uHdPO8e5aAjVlMGY+CxJ7NUhR2GQbxJ0en/H4KBfmjtnPueLu08GgfRYdukUVj+ggb4ezmQYFt7+zjWBfEzOG9WLOqFhGJYQ0+yGpzzJn8diKvdRYbbwx4hC+Byph0l2AnvZRP6Vgw5ECrl+8kUc/3cs/5o3QV0jPeRXemwtf/h5mdWyx1/LULJ5ZdYDsokpign24/rwEDuWWsmLXCVTteG6b1LfubzYlUX9e/ecsmJbUbBrEgmlJDf5+QA+0H7xkIGXVVj5JzeT3S3fxx0/3cPHgKOaMimXigIgWL/hsNo3vDpxk8Q9pbDhagJ+nkRvG9ebXE/qwLb2wyfm8PFy3RsGu/s8tzN+TwvIaBsUE8sbNY/Brpu5tezV30XaimbQZ8g7Cj8/CsHlw1WsA5GZl4f3qWDy/ehCS1jm/0OpM1vg9BeDnL2DXEnKKQ0lq5/twc9r8e2rGRYOiiAny5u0Nx5g2JLrlncP66e/VZSchIKrlfdsgLb+s1fxmuzF9QvnuwAFOldcQ6ufZ+hNakbrmY6YZdpN33iJ8fEPbd5CgOMo9gulV8TOaprUrXbGJ4AT9q28YVBTo/z9kToPPp/c2pZNfVsM9Fw2o29Y7TJ+1P1ZQzoj4lu84tFl2KkQNBg8viivNlFRZmk+Z8fLHOOgyZhz5jlOX/5m/rz7KJc9+j8Gg6u5uuDIFsatJ4NxJ1h/Op7TawmXDTqc3jEoIoW+EHx9tPd7hwHl3VjGDogPaNwvM6bJ0Zd7RBFQ5aJIQFMex/HK+2pvDHZP7te8DXNNgy+sQPRxiRzM7Tn/DsQcJ+0yDsGgGEkpTgVYK9LdFS7PoXRw4a5pGmJ9nk8DZx2Tkr7OHEOLvxfLULJZuz+S9TRkkhPoyOzmWK5Nj2Xm8qO5nFxXoTaifiX0nSknpHcIzc4fR54M/QmwKxI12eO5x/cL43dQBvLDmEOcPCOPK5DgYcDGcfz/8+E/ofT6MuKZd/y5HedvPrDqAp1Hxq/GJ3DKxD72Cms5eNA7uW2Lfr7lA++4L+5N6vIjlqVms2JnN57tOEObnycwRMcwZFUv19g+I3/53IrU8clUE38b+hteLUziaX06vIG/+cNk5XDMmoW6Rlf3Wsv18KEhy8cKcxj+3/LIaFHBtvXG4Ukywj8O0Foez6JoGK+8HT1+Y9te6zXGxsTzp/SseKn4etr8FKb9y+Ti7HUfvKZZKtDWPkV/2tMsrarTp99QCD6OBX4ztzTOrDnAot5QBLQX49UvSuShw1jSNtLxyUno7F7TWz3NuNdBvTr01Oxdi4BRBREx10HzLWUpxKmgIg/OOUFBeQ7i/C+p1H1qlz+7evUPvZ/Dm5XDgS6gsBJ8QqsxWXvn+KOP6htX9TOB0iUKXzzjXLgxkyJWAPtsMEOdgxrnOsKtRe5ZyY2QaMx6cwsSnv6tby2XXmWtC3OksmBroHr7cnUOAtwfj+4fVbVNKcfXoeLYcK+Soo5bUTrLZNPZklbQrTcMuNtiHAZH+vON7o36LqD4P/ZbRaz8cxWQwcPP4xPadJGMjnNwLY+ZD7VV6/U5n3z88i8PGfhTs/Y6Sdt5ed6gbVwH4aGsme7JLmDWiV4NWyU/OGcac0fFckBTJ89cms/WRi/n71SNICPXlX98e4oK/r+W+JTvqKjzklFSx70QpV46M4cPbx9GneJO+WOm821s8/91T+3NuYiiPfLKHtNo2uFzwMCSMh8/vhbz25Vk3d2s5zN+LRy4f7DBobo/GnfLqvyErpRiVEMLjVwxl0x8u4rUbUzivbyj/25TB4pf+xpBtfySaPAwKepHHlZl/Y3L1d7xwXTLf//4CbpvUr0mwWv98d07px+7skg69dhtz9HPTgP+sdU+d7QXTklpMj2lg14d6SsJFi8A/ssFDNUPmsUkbhLZ6kd49r6dr9j0lC5vm+ooabfo9teLaMfF4ehh4e0N6yzvWBc6u+9vLK62mvMba6sJAu+FxQXh6GNpfz7nRmh0PrASpCti3vH3Hq2WJGskAlUlGbkGHjlPn4Nd6dQ/vQP2zcfpTUFUE3z0JwPubM8grrW4w2wzgbTLSK8jb9YFzYRpUFdelpWQW6hdtLeal97sQvINhz8eE+XtRWdP0/R96RodRCZw7gdlq4+t9uVw0KKpJXtdVo2IxGhQfbWt/EGdfGDg8rv2BM+jVNZ7LTaZ6+nN6C1L7AqCk6eT3vYKPt2UyZ1Qske39UNiyGLyCmq3a4ONpJGTwBQyyHeSJT1Pbdw4HKnwcz1Q0t72zpBeU89iKvYzrG8Zz1yQ3GwCCnt87d3Qc784/jw0LLyTQ28Nh6uHmY4V6tYVNr4J/FAye3eIYPIwGnrt2JCYPA797fzvVFqu+on7u63p77o9uhpq2vyk39+boqsWwbeXpoeeCvvSL0Wx55CIWei5pUovYV9Vwu/k9Zo2IcerOzU3jEzEZDSxuYzORlnR2O+vZybE8OWdYk4u2JjNCFadg1cMQNwZG3dzkOJOTIvlDza/Rasrhmz+6ZazdSjOLucx+elUjVy8OdPr35IQwfy9mDo9h6fbMlicoghNAGV1aWeNo7cV5SzWc6/PyMDIyLpgt6e1cIOjgzoBRMzss89YWXr1H46FslB7b3qHjAFCUAXn7YcC009uih8Hom2HLYqqz9/DyuiOc1yeUsX3Dmjw9IdSX9AIXpzfWWxgIp5vutDjj7OEJg6+A/Z9DTUXndjjsZBI4d4KNRwsorjQzfWjTQC0y0JspAyNYui2z3S08O7IwsL4pSZHUWG2s970A7tsDi4qg31Q49iPvrT9ItcXG/InONTxpouwk7PsURl4Pns2/aUYNvwgvZeHYznV8s8817V6fNl+DWWt4wVKlmXja3L40BFewWG3c9+EODAbFP+aNwLDnI6cXL0YHeVNa5XhBaXZRpV579dAqGP0r/c2sFTHBPjx91XD2ZJXwty9rZ5gDY/R855P74YsFbfq3ma02vLpoIZ0zgnxMRGv5Dh+LbGa7w30DvLlqVCwfb8skv8w1zVC74sOmpVn7OqsX6beNL3/WYQ7z2L5hZBrjWR95vd7ePe0Ht423W7jwUX0hbX0mH/YPuQ9w/YwzOPl7ctLN4xOpqLGytKUJG6MJQnq7NHC239Wqy3F2YtH2mD4h7M0qpqKm7YvotWbuDDS33VmhA/Sa+LYsF0zwHFylfx04reH2Cx4BL38Klj5AbkkV91w4oOlz0S9C0p1ow94m2alg9IKIQYA+4+zv5dF6utiwq/VW7Qe/culdku5GAudO8MXuHHw9jc0Wnr86JZ6TpdV8f6h9tzg7ujDQbkyfEHxMxobtt8+/D8pPUrbxLS4aFEX/SOe7EjaQ+g7YzDDmlpb3SxiLpgzMDDzCQ8t2UeCCgOStsjEUan5UaSZsGlg1xVEtmrfKzu3wsdvrP2uPsD2jiL/MHkpMxgqnu0bZtRhgbX5Nz5drQ67pJUOiuXl8Im+sT2O1/YKl/4Uw6UHY8S7s+J9TxzFbbdzzQSpVZhumRt1+us2bZlUxVuU4R79Itdw+vrH5E/tittp4+6djLhgYPHjJQBqvNeryn1vGRj13eewd+kyYA94mI2P7hvFE2Qx9pnLlA2BxQ4ms7mL4PIg4B5QBUPodupkvsDPkYgCiglyQ9+pGw+KCSE4I5u0N6dhaKk0X2telTVDS8svx9DDo71NOdssbkxiKxaaRmlHU5vPl4rhBV3PbneUdGkc+wfgXuKBb3qGvIaQPhPVvuN0vDPPkh4gp2Mid0T8zrl/T2WbQm4rllVa7pDpXnewdED20buIls7CSuBCf1hdC9h4PAb1g98cN7pKAnoHy+BVDzvj8ZpDA2e2sNo2v9+Yw9ZzIZstITT0nkjA/T5Zsad9VcEcXBtp5eRj1snQH6wXOiRPJCxrGL62fcvvEhPYd2GaFrf+FPpMg3PFVcx3vIFT0MGaHHqO40swjy/d0qCakxWpjuGc2kYYSHrfcSN/q//F3yzUMNhzngkA3FY5vxY7jRTy35hCzRsRwxcjYNnWNsmvuav6hqXGw4z0YMhsC2paKsnD6OQzuFciCj3dyorh2PFMegsSJ8Ond8PekFmeGLFYb936wgy925/DIjEE8M3eES24tu1RNBfzvWoxYqdYaBs82TRFMKexzvrFpvwh/LhoUxdsb09s1I9ZYsJ8nmgbBPqbu8XOzmuHz+yAwTv9baMHkgRHsz7eQN/GvkH8ANvTgFtzFWZC7R59YWFSk36EbPo+c4iqMBkW4X/cOnEGfdU7LL+eHwy3cZQnt59KSdGn55SSG+erpZE6+743qHYJSsLkdec4v1sxoMvQKzZMnazq++DzdK4nosv0dO0hNBaR9r882OwhKP9Qu4YAtjrstb6IsjieR7GkvGa6adbbZ9MC5Nk0D9MWBTjWhMRhh6FX6xUBlYd1dkv/deh6aBgZXVCDpBiRwdrPNaacoKK9pUE2jMU8PA1cmx7J6f26bZ1htNo29HVwYWN+UpAjSCyrqbqlZbBrPVs6gt+EkKeXft++gh77WZxPGzHdu/97n43cylQcuTOTLPTnt7oxUWWPl9ne2cYFtAzZN8Y1Vry7xrvUiSjUf/hz+dbuO2xEVNRbu+3AHUQFe/PkKvT5mexYvNpfzeDnroLoEzvtNm8fmbTLyr+uTqbbYuOeDHXoZIYNRz1vTzFCWQ3MzQxarjXs/3MHK3Sd4+LJBzJ/Y16W3ll3CUgNLboSMDRjmvs6u0U+QQwQ2TZFDBDtH/glDXAp8dBPseN/pw/5mcl+KKsws2eKgm14baJrGs98cJC7Eh80PX9Q9fm4bX4KT++Cyp8Gr5btN9so8q8wj9FbG656GwmOdMMgukPouaDYYdWODzTklVUQGeLW5FXpXmD60F+H+XrzV0t2S0L5QUwrlrlnwmZZffjpNw8n3vUBvE4OiA9tez1nTmOG9CzMGcrRgbJoi0xbOQvN8tgZe3I7RN1QQNJhY63Go7sDi4GM/gqUKBlzS5KEai42X1h3jf6F34l12HDa86PAQ9pJ0LstzPnVE/53XBs6aptXNODtl6FX63eX9K+o2jesbRr8IP97e2MqC1DOEBM5u9uWeE3ibDExJaqaZQK2rU+Kx2DSW72hbkJh+qoLSdnYMdGTyQH21/LoDJwH4am8O75cMpSygH+rH59o387BlsX77Juky5/ZPnACWKm7tU8iohGD+uHxPmxeVFVXU8IvFG/n2wEluCt7NqfBReAbrFy+l+JLW93pis7/R69J2or+u3M+xgnL+Pm8EQb4m/efZXM63T8slm5oEpiN6weZXIWYUxKW0a3z9Ivz58xVD2Zx2in99W9tBcP3zTXesNzNksdq4b8lOPt91gj9cdg63TmpnHrw72azwye16a+SZz8HQqxgz63aiFx3G8FgR0YsOk3zlfXDDJ/oM+/Lf6CkvThjdO5TRvUNY/GNau9cpAHyzL5ddmcXcfeGADncGdImiDFj7FCTNgHNmtLp733A/4kN99FSvS5/SF5Z9scDlDTS6nM0K29+GvhdASGKDh3JLqtyS3+wOnh4Grj83nu8OnGw+6Kpfkq6DrDaN9IJy+oTXXoA11y0voOkk07l9QknNKMLcltfXlsWMs23nL5YbGVv9En2r3+P8mhf4xjjZJalPNZEjMKBRebwDec6HVumLsHtPaPLQx9syyS6uYupl8+Ccy+GHf0JJ0/ggoa6Ws4tmnBstDCyuNFNWbXE+cI5J1u9U7P6obpNSihvG9mbn8SJ2ZbY95aa76Qbvzj2Xzabx1Z4cpgyMxNez5brHSdEBjIgP5qOtx9uUmuCqhYF2CWG+9A33Y+3BPDRN49Xvj5IYHoDv1AcgdzccXt22A546qj9n9M36YhOnBjEOUBgzfuIf80Zitmr839JdTv9csosqmfvyBvZklfDmzFBCyg4RPuZq1i+cytf3TQLg596/AA9vx0Ghm6zZn8t7mzK4dWJfxvcL12+JrXxAbzJgaPT3oQxQWQArH4RmbtE1cfQ7yD/Yrtnm+q4aHcec5FheWHOIjUcLWpwZslhtPPDRTlbszGbh9HO4bVLbW9i6nb3+8N5lcPHj+t9ic7z84folerD4xYPwwz+cOsVtk/qSWVjJl3ty2jVEm03j2dWHSAzzZU5Xz8yD/jOzLwqd/jennqKUYvLACH46kk+NXwxc8Af9blO9mace4fAaKMl0+Hfkjnbb7vSLsb0xKsU7zZWms7ekdkGec1ZhJWarRl/7jPOFjzYtfQr6e19VcYNNYxJDqTRb2dtMp8sm8g7A14+QHjKet60XE+7v6fLUJ8/4UQCUHNncvgNoml6Gru8UMDX8m6mx2Pj3d4cZGR/MpAHhcMlfwGaBb/7U5DCB3ibC/DxdV5Iua7tegjZcv7g4fsqJUnT1KaUvEkz7AUpPvx/OGR2Hr6eRd3vArLMEzm60PaOQk6XVTB/mXK7pvJQ4fs4prQuGnbEnqxhPY8cXBtY3aWAEG44UsO5gHrsyi5k/sQ+GYVfreY4/Ptu2g239rz7z1OiWZot8Q/U2n8d+oE+4Hw9ddg7rDubxv80ZrT71YG4pV/3nJ3KLq3j7lnOZbN2gP3DO5QAMjApgQKQ/Sw/U6GPa9QEUdewWuzPyy6r5v6W7OCc6gAcuGajPWq24G7a+DhPugdn/OV0CMCgerngJxt0FW16DNy7VZ/9as/lVvQ3ykJZL0Dnj8dlD6R3mx70f7MAa6PhDRguI5sGPdvLpjmx+f2kSv5nsgqDZ1a3RNQ2+eRS2vQkTH9B/1q0xecO8t/TueGse1z+sWrlou3hQFH3D/Xjl+yPtysn/am8O+0+UcM9FA/Do4FoFl/j5czj4lZ7XHLdsPqUAACAASURBVOx8c6bJAyOpqLGyNf2UfgEXNRS+WgjVpW4cbCfb/pb+OnNwBy23uMrlzU/cKSrQm2lDo1my9bjjHH0XlqQ7mq+nNPSx13AePg+mP316h6B4mHAvlOXC+9eD+fRdxjF9QgCcq+dsqYFlt2Iz+fHL/Ju4dEgvtj5ysctTn6Jje5OthWLNbGdJuryfoTjDYZrGJ6mZZBVVcs+FA/QFeaF9YPxdsHsJZGxqsn9CmAtL0mWnQq/hellSnGx+0tiwuYAGe5bVbQr0NjE7OZZPd2RTXOHCPg1doBu8Q/dcX+7JwdNoYOo5ka3vDMwcEYOXh4ElW50P5HZnFnNOrwCX3tr19FBUW2zc/N8tGBR4Gg366trxd0H6eocvXIfMlXo1jXNm6OXN2iJxAhzfDJYafnleb87vH85fV+5v8c1hW/oprn55Axabxoe3j9NrXu5foacu1Pvwv2xYLzYfO0X+8Nv0Dc3kjrmKpmksXLqLkioLz1+bjJfSYPmd+s9m0u/hosf0DxF7CcD79sDI6/TubPPegYLD8PJEfXaiOaeO6mWNUn4NHh1fmOTv5cG/rkvmVHkNr3vegNZoZkgDqsuKyNi5lgXTkrhzSn/HB2oLJ1fZt8mP/4SfXtDz66e2ob6w0QRXvqL/PNc/p89Y25q/TWwwKG6d1Jc9WSVsONK2pghWm57b3D/Sn1kjOvih7ooLj+pS+PL/9KB37B1teur4fmGYjIp1B/L0D97Ln4WSLD3loycoOaF3dBt5fZNSj+XVFkqrLUQGdv+FgfXdPD6RkioLy1MdpAkaTXrw7IImKGmOajhHD9e/zntbf9+7+DG48mVI/xGW3qJPMKCXfkwM82WzM3nOa5+AEzv5b9j95NqC+MNlgzo8dkd6h/qx29YX3/ZW1rCXoWsUOJutNl787jDD44Iapnief7+exvLl75u8FyWG+blmxtlqgZxdjRYG6jPOTi0OtAsfAL1GwJ6PG2z+5Xm9qbbY+Gib+yer3EkCZzfRND1NY+KAcAK8nUtRCPQ2cdmwXny6I5sqB13XGrPZNPZkF7ssTQP0lr/1O0rZNHj0070sT83SZ2h9Qp2fdd67XK/96uyiwPp6TwBLJWSnYjAonp47HKNSLPhol75orZHV+3K5/rVNhPp5suyO8QyOCdRnkrO3w+BZDfadMbyXfif6uAcMvxa2veXWbmfvbz7O6v0n+b9LzyEpwhuWzddnuqc+AlMfdriaus7gWXDbWn025n9Xw5o/132YNLB5sb6Qb7Tr2h0PjQ1i4fRzeCJzGOsG/rFuRlwLiuezyDvIsQayxPsJfhvpgpJM0K7qIi3aslh/7rB5MP2Zln/OjhgMMOOf+iz11jf0vGdr85UzrkyOJdzfi1e+b9vs3Oe7sjl0sox7LxqgVxtoL1ddeHz3pJ5LefmzzqdX1fLz8mBMYujpyjzx58Kom2DjfyDHRX8nXWnHu6BZ9X9TI7kl+gzpmZSqAZDSO4RBvQJ5e8Mxx3dLQvu6ZMY5Lb+cAC8Pwv3rXXDk166jCKtXbWnYXLj0b/pdj5X3193tSUkMZeuxUy2Xzzu2Hn58jvyk6/jz4T7cMrFPXQ6wqwX5mjho7E9wRXqT1BKnHPoaooZBUMOL5U9Sszh+qt5ss52Xv55qdmKHXjmpnoRQX7KLK/UGVh2RfxDMFQ0C5+OFFQR6O1HDubGhcyFrW4M0n8ExgaT0DuHdja2UQezmJHB2k12ZxWQVVTK9hWoajlydEvf/7Z15eFTl2f8/T/adkEnClpAESFhFVKAqi1WgVrto69ZFq75t7ebb2lZbtf7U2vpqpZvWLq/dtH2rLVpFrSAKLuCGAkFAgQRIIBuQPWSfZJ7fH89MkpnMJLNlZpLcn+vKNZOTM+ecmZM5537u53t/b0519rDpg+G1kkcb2jnVGbzCQDAtfzutzqNZR3954pLN9GvJRjjx4fAbe+9P5oJYsNL3A3EUSxx9AzD+xHd9ej7vljfwF5dObet2VPC1/9vJnMmpPPX1c/q1WA5t5VznwNkh13hhTw0sv8lUNW//ve/H6AVldW385D8fsnxWJtd/ZKrpxPfBM7DmJ7DSy8YilpnwlZfhjGtg28/h75eahjIOulpN9nrepZDm2//bcFy/LJ9Vc7L5r135LGn/NQWd/2B+48/5zrEVvHTuP4jJOcO8pzd+FXgRWDBbo+9ZZ/ThRRfBpb9z27TDK5QyN6tVd5qW039cBb+c7zajmxAbzXXn5vF6SS37a7zTYvb02nhwcylzJqdy8YIAz93muwIfeNS8b74LZ11ngl4/OK8oiwPHT/UX9K6+29QT/PGC4ElwwoHNZooCC1b2a38HcHyUBs5KKa47N48Dx0+x3Z0UwjIT6o8E/P0uq2ujICvZORisLzWaZkcRooOzv24yrDsfhddM2+ml+Rk0tls57KnFfWczPPM1dEYB366/nMyUeL51fhBmwoagLm2eeVLzvm8v7Gg0/uhF/dnm9cVVnHv/Fn7w1B5ioxUtHW7kDKddAbkfgS0/dgrW8zOT0Lo/O+w3fYWBZ/YtMo4afgw+FlxmHh/5qNP3/ppz8iivb+eNoWwQIxwJnEeIDftqiIlSrJk7yafXnV1gITcj0Su5hkMLHczAediWv0u/CrHJwxfVVe+Gqh0m2+yPd2OyxXQtKn+zb9FlZ05jzbxJ3L9xP0vv3UzBrS+w8O5N/OCpPZw708LjXz0bS8qAadL9z0P2fLc3OYdc42Rcrsnqvvsn/7IGQ2C1W7TFxUTx88/MJmrdNSaL8vGfwbJv+7ax2ES45GG45LdGwvKHFfDKT83F6L5pxoLO1UA/CCilOH9ONjYNtae60EB7dy8xUYrsSVPhS8/B/M+aznLPf8f4/vpKW51p54yHG7OKMoMwb4skD26EZ74O+cvhikd9zpq6ZcX3YeHn4fhuUxjmIaN79dl5JMVF80cvs87rd1dzpK6Nm1YX+Wdh1nQM3v4d/OUitxX3gP04b4IDG4a2zrL1Gs/mJAusHlyE5C0fnW135imxD+4ObTb2VL3dBE2CEw6OvGo+bw/FpY6M86RRpHF28OnTpzEhMZa/vV0++I99lnSBBTpOVnQO6kogPW9QcRxgBqtnXA2v/wze/SNLCozLkEe5xgs3Q0s12xbcy1uVXfzgwtmkxA9dlB8onVl2qUn1bt9eePhVM3Nhl2msL67itqf3Ut1k/oesvZrbn9lnZnoHopQp1m2rM3aPdqZnmM81YJ1zdTHEpTjdS4yHsx+dS4++aa7dXS0M/N5fzDYsyXH8fRQXCUrgPAJordm49zjnzso0lmM+EBWluOKsXN48VN/XH94TI1EYOGzL36QMc+PY+yQ0DvGPv+PPxmbn9M/5fzD5y6Fie18wppRiZVEmvRpO2oO4ls4eopXiktOnkjzwInnqBBx7e5BMw4FDrvHiB8dNZqOrGd77s//HOoD1xVUsu/8VCn+0kfcrmrhiYQaTN1xvrIc++SuTTfGXM66Gr2w2F92ta+3T8nbeenBEgpHfvzZY39hj02YWIjYBLvuzCSx3PQaPXwmdXla+dzbDK/fCg6cbv+Dpy0xmciDR8cby64Xvw2/OMhmooYLzsm2w7lqjr/v8E+5vyP5in/1wwtphNIdHXofmKtITY7lycS7PvV/dP9j0oDu29tp4aEsp86emceF8HwbYtQdh68/hf8+DX58Gm24zN6d4D10PYxLN9/Wfn4cHCuDvnzHSCccUquP47skwU6tzPw2JE334YJwpmpTC5LSEfrnGlnvsQfMAApHghIudjxqpmr3Q2JXjzWZgN1rs6AaSGBfN55bksumDE4OTJxn2xEMAco1Oay9VTR1uAudDnptiKQWffNDMGm24hfwTL5OZEs+O8sbB6+59Cvauw7r8Fn64PZ4F09K4/CwPdndBJDN7ChU6C1uVjwWCpS+Z71jOEsDM9Ha4yDP7ZnpdmXqG6dD39sN915SiExsBKK8LUOdcXQxTFvXN0GmtqWjo8N5RYyBb7jFe5wOxdhD76k/53NJctuw/QZWHRF2kI4HzCPBhTQvHGtq5eIFvndscXHZWDkoZH8ehGInCQK/6y5/zLTOS9FRU19EEe540WrXEdP8PJn+ZsWobMA32h9cGX7x7tbHycuLAfwBtGjG4wUmuMXURzFxlgjfXqW4fcWQOHBeEJDpZvfvb6MOvmmzx4v8KaPuAaXscHTd4+QgFI8POQkRFmezQp39jumD95eNDyyu62+GNX8OvF8LWB2DWavjmdvivDWYbTu4iD8N/74Srn4aUSSar/ZszYdffTQA9MCj9eSH832WmAv3qf0N88AaUgOf31NEIf/s0/Goe/M9Ubj/2VR6MfpAj/7rVZHCfvdGt7vipnZUca2jne2uKnKev3QXa1cXm3D68BH67FF75icmkr7kH/nsXfONN+MQvBtt7xSbCpx+CH5SZ2YGlN5j38eKt5nNcWwjrv+E8AHv/8YAGYA5bum2ldcbXOpgSnHDRehIObrAXBbov/jvR0klKfMyIZzlHiqvPzsOmNY9vd3Hw6fNy9r9A8FhDO1rjHDjbbEaqkVnk+YXRMXD5XyD3I6inv8rVk8oGdxBsqoD/fA9ylvI72yXUNHdy5yfnh6QJTV5GMntsBfRW+eDlbLNB6cvmuhdl7rXDXmMHsmedGeACjmtKysvf46r4twPrHthrNXUIUxf1LWpo66bD2utfxnmI7/3nl5ouxI9vH51Z59H5DY9wNu49TpSCNfN8k2k4mJaeyPJZmTy1s5LvrCp0ewHQ2hQGfup0H90qhsFh1bN200GqmzqYmp7ILRfOdrbwmTANTr/K6P1W/gBSXJq7vP+EKexb/OXADsahcy5/o6+hh9cXmP3Pm0xJ9jyPm7/4tCk89EopJ1s6yV7xfXj0YtMRbOlX/T7ktZsOsqb3dX4Qt46pqg4rMcTSwz2xN3HXGVf7vd1BeJyWD34wMjU90W1mYNDsxJlfMsHuui8ZLfCSL5v/keZK0+zg/Nuhu81kyltPmGnKC+4w2WEHC680P67MWgUzLzDT/q/eC8/daOQhXS392czWk4Ayg5OkoZvH+MWEHOcA00HqFOPAUX8I6g8RV3+IJc0fkFH9LtS4ceKwdqCf+28stiX8Nn0CF1TvgqYMSEiHEx8YWUqvXZbSXAFP3wBoYwuWv8wEv+6cahyf25Z7+j/zVXf2L59xnvm58F7TRvnQZnjpDuMP63J8bLnH/Xnwko/OzuJfOyoormhiiafPLTrOTDknZ/q9n5Cx+x/mcxrCA9w0PxldjhoDyc1IYtWcSTzx7jFuvGAWCY4EShAs6Y7UGgnBDEfzEzD/Ez2dnjPODuKS4Av/hL9cxDeP38XL7T+iuukcc/2x2czAT/dyYvVD/P7P5XzitCksLRiB778bpluSeM02g080v2sG0N7M1FTvgvY6KLywb5HX11gw380e54ZgytrB96L/yQ/r3c+wesXJ/ea6E6ijhgNP3/sJOeRMTOKCOZP413sVfHtVIfEx0YPXi2Ak4xxktNZs2FfD2TMsznpbH7licS5VTR285cHa6mh98AsDHXjVKnnZTUZzuv0Pzsu1Njf+nCVOI1e/SMk22Yij/TrnYaUkAO0NUL7NyDSG0Fc7yTXyzjVFF28+6J9O187ilpe5P/ZP5ETVEaUgXvVgJYaGtu7hX+wLnrpueVoeAF7NQjiYeT58+SUTzL56r3Omdf03TVORjJlw/YvwxSedg+bhUAoK18BXX4XP/9NIPVwlAGh46zc+v0evcNewITbRZH1nnGcGCh+/D774JCeve5u5XY+i8fD/19NJUU8JH7O9gdr2C9h0Ozz7TXjnt/1Bcx/a3JBvLoVrnzcDO0/2jq62hp6C34wCsx1PuvEAB2DnzsokOspuS+fuc4uONYHo/66EivcC2teIY7MZ5528ZUMGecdHUddAT8zKTqa+rZs5/+9Flt3/itHYxsQZO88AmqD0WdFlDgjA6t04angicSJc8zQkpvNo3M+o3/xLu7xoornWL7iM+97pxKbh1ovm+H2cvpJnSWKvLjC/eKtzLtlkZmxnrepbdPPHBmfdPV5jPXw3s211HAvEks6lYyAYRw2A3Aw/Ms5uG9woOOe/AbjmnDzqWrt50c+mUeEkKIGzUurjSqmDSqlDSqlb3fw9Xin1L/vftyul8oOx30ik9GQrR2rbfHbTcOVj8yaRlhDjsUhwJAoDfSKz0Mgg3vujs6a17HWTefPHgs4dectM9bHdgs2rIO7gRnNTnjv06NtJrqGU0ek2Vxi9nJ/cGreOJOUczMWrHm6Le9LDK/zEUxC36s7g7gczkLrvs6cxLT3Ru+5b2XMHa5UB0Ca7eP0GyDvH/wNSCmZfNDhT6mCkJAALr4RPPeQsJfnUQ26D0wXTJrB01iSOY3G7qeNkcsvUR4m5rRzubIAfHoVv7zbbdUdHkymYDTYjNACbkBjLmdPTjc7Z3ed2ye/ghldNx8y/XgTbH4nc1tzlW6GxbOiOk9ibn4ziwHl9cRWPvVXe93tVUwe3Pb3XBM8ZMwPKOJfVtZKZEu9szeqwohtKqjGQtKlEf+kZEulkwT7n+o7e99dh27OOG1bM8E+P6yeTUhM4GGXXgFd7Kdco3QQ5S51mxQrtdUrpSbHDX2M9fDdb4iZR0dhu5FH+UF0M8ROcHE4cGedpHhJWQ+L6vU/ONrNM238HTcdYMSuTfEvSqOwkGHDgrJSKBn4LXATMAz6vlHKdH/8y0Ki1ngX8CvCuh+soZOPe4yiFb8U+bkiIjebSM6bx4gfH3XbZGYnCQJ9ZfpPJ+u18tH/Ze382BTTzAu9eB5gCwa4WY8qOl0Hc/ufNl3XAyNkTfe4aLZ1GOjBpgbFWG6LZhVt6uuh57y9Mxn3l+SQPy/3GhyAuGHg1CzGQUzXul7fV++ey4o4QZt378DajC3xt5Uzu676SnmjnYMoalcB93VfyXYe2OSrK1AJkFIT+PY3gAOy8oiz2VjVTe6rL/ec25XT42usm87bxFvj3V4Z2/fCFYHag3PmYkdEMMRC32TQnT3WNSkcNB6ZAzYMVqcPL2c/BTVldW3+rbQd1JeZz9UGqE509m57oJJSLA090bwe3xa3jGx8NQvdSH4iKUqRnZHMyZqp3gfOp46Zmp8i56cmW/SdRCjZ/77zhr7EeWpUfLLoBa6+mprnTzYu8oLrYzBIPuD5XNraTnhTrdS+KQQz83t9SCtdvNJKWv15MVFMZV5+dx3vljV7bd0YKwcg4LwUOaa2PaK27gX8Cl7iscwnwmP35U8AqpYJ194wsNu6rYUleBtmpgV9Ar1ycS3ePjeferxr0t71VzcyeHNzCQJ+ZdhYUnAdv/9ZM+bZUw4EXjPNDsNwMBuqc7QwZxHWdgsOvmGy4F/9iTnINpWD5d6HuIBx8wbvj62o10oAHTyfmhe9i1e7LBtRIBD4+BHEhJxQBYAiz7v6wojCT0kkX80DsN9H2AY4tLYe79Q3U5l/CuTPdBAyhfk8jOABz2NJtKx2iuVDiRPjcE6ar4wdPw59WQW1JYDsOZgfKtjozEF/0hSGvafVt3fTY9KjOOA9ZP5IxwyQw2n3riumgrK7djaNGqZm59DEUmGBrcrt8MvXOzkohIi8jiQ/VTO+kGqUvm0eXboFbDpzgjNx0Mr2Rd7p+Z1NMki6nx2Tg/eog2NNl6itckk0VDR3k+qNv9kTOWaZIubsN/voJrszvJD4matRZ0wUj6poGDNQTVNqXuV1Ha90DNMPgOUyl1A1KqR1KqR21tSPXyW2kOFLbyoHjp/i4n24arsyfmsbcKWk86eKuobVmX1VwOwb6zYrvQetxUxC48zFjP7M4eN3rSJtipgkH+DkPSckmoxEdRqbhwEmuATD/M+Ymse0XQ2dX2htMh7VfL4CX7qBjwkyu7f0Rf8++JaKDuZARigAwxFl3X1FKccPKAh5pWsyrF22Bu5v441nP8Y+Os/meG00jEJ73NEIDsHlT0shMieu3pfNEVBSsvBmuecYEqn88H/Y97f+Ot/w4eB0odz9uPKjddAocSJ+H8yguDhyyfsThhe+Hzrml00pdaxcFWe4CZy9lGgPoTvYggxzJmaYhmG5J4r2uPGg+NrzXdekmSJtmZjbtnGjpZE9lM6t86fkw8Dt7cwmceS1TSv7ODFVNuT9ezic+MP/nLoGz3x7OQzF1EVz3H7BZSfvnJXx1dhfri6to6fS/tijUBCNwdjdcdI04vFkHrfUjWuvFWuvFWVlZbl4S2Wy0i9yDFTgrpbhycQ57KpudpjKONbTTMkKFgT5TcJ4xsH/h+/D6/caqqXJHcPeRvwyOveW+1bQr+583Wiofup45yTWiok2L5epiOPLa4JVbakyzjl8tMO93+rnY/msz1/Xewa7o0/nUNTdFdDAXMkIVAEZy1h345MKpTJ2QwB9eP0JrVw9/eP0wKwozWZI/RNV/hL8nb4mKUqwszGJrSS293rTXnfFR+NpW44Tz1PXw4m3eF+p2t8GHzxkHkmDZ32ltvMlzz4bsoQvO+gPn0Ztxdlc/EhutTP1InyWd7zrncnthoFPGubPFJFz8aNqkVt9Fu3a24+yNTkSFKTmRl5HEzp5888tQWeeebjj8milwHpBl37LfNAry14ULMO5EMYncEfuEf5Z0bgoDtdb2roFBDpwBJs2H614AFcVNVd8lz3qEZ3YNnlmPVIIROFcCuQN+zwFcvbL61lFKxQATAA/tf0YvG/fVsCg33ePI3R8uXTSNuOgontzRf9HfUxnmwsCB7H3S6FkdhVo9ncHvCpa33GipT3ww9HrWDjMVNveTff6Y3uAk1wA4/fOmSOLxq/o1km8/bDyEH1xomkfM/SR84234/OP8s2YS28sa+NHFc8lOSxgzgU/AyOdAbHQUSwsyeLesgQV3baKx3crifP+bi4w2zpudRWO7ta+YeVgmTDM31I983fiqP/pJ2P6/7vXK7Q2w+wl44gvwwExYd41pLBHrYWo5Jt40RvKWo2+aQudhigJhQLvtUaxxdq0fiVaKGZnJRgqXnmecIPwInB2OGk4a53ofCwMHEHfG5/hF/Leo0pnYtKJKZ7Lr9B+H7fqSZ0lmny3f/DKUzvnYW6YD4wAbOoDN+0+Qm5FIYXaKhxd6QUo2auX3uSBqJ4nHtvr++upiU5uUPr1vUW1rF109tpErtsyaDddvICYukXUJ9/LOm1vQkVog7EIwAuf3gEKlVIFSKg74HPCcyzrPAY65rsuBV/Ro+YS8pKKhnX1VLVx8WnCyzQ4mJsexZt4knimupLvHFG7sq2omNlpRNDmAL1qwCEVXsHy7zvnoMHKNQ1vA2ua1TMPBILnGh8+Ctd1uC2bXSG76kWm6ccbVpiHHZx+BSfM40dLJfRv2c84MC1ctyR1yP8L4Y31xVf+AzM4fXjsyuJXuGGVFYRZKYWzpvCUmzrQVvuzPULXLdGZ0tTV8+COwdhas/zrU7IYzrzFWfTcfgk89OFgmFBULvT3w+3Ng/3+8O46dj0LCBJg/fKHzieZOlIKsACxII4GB9SM3XjCLkpOtxl84Js7MGvnRBOVIbRtKGUlDH746agxgfXEVj7UuZVnXQ8zo+gfLuh7iS+/lhe07Nd2SRCtJnEouGDpwLnnJdEKdcV7fovbuHt48VMfquZMIuOzrI9+gNmYKl578rflf94Xq3Sbb7FQY6PBwHoGMswPLTLh+A1EJafys9Q4qnrkzeEW9I0jAgbNds3wjsAnYD6zTWn+glLpHKeWIYP4MWJRSh4DvAYMs60Y7G/eZoOuiBYHZ0LnjisU5NLZb2bLfZEv2VjUzZ3JaZJiGh6Ir2IQck/EYUCDolv3Pmyrt/OU+78JJrrHlHqP3ciVlkmmZnWF8O7XW/L/1++jutXHfZ08L/MInjDnWbjpIpyengnFARnIcC3PSea3kpO8vPu1y941sbFZoOGQkVV99Bb77AVy8FgpWmk5z7mRCl/4OvvmW0Zf+64umm+NQDh7tDUb6sfAqtw4Grhxv6SQzJZ6Y6LHTGuHys3LQGv7tqLGx+GdJV1bXRs7EROf7VV2JaaoyMd/n7a3ddJAeF+lPOL9TORMTUQoqk2YPHTiXbjL3prj+zPsbpXV09dhY7Yu+2ROxCbw+/UYKbEfRu/7m/eusHXDyQzf65gCan/jCxHyiv/wi3SqG3D0PBaeod4QJyrdca71Ba12ktZ6ptb7XvuxOrfVz9uedWusrtNaztNZLtdb+G0JGGOuLq1h2/yv8z4YDxEYrdh5tDPo+VhRmMTktgXU7KiKrMBBCZ5+Vv8JknD3ZxPV0G//mOZ8wzRV8xEmu4Snod7FY27jvOC99eILvriki37ViXBDwsZXuGOWjRVm8X9FEoz9NgFo9SCtsvbD6LuPs427A6k4mlDUbvrIFln/PdAj9w3KoeNf99t//p5lx8kKmAXCipWtUO2q4IzcjiXNnWnhyZwU2mzY653rfLenK6trIt7gpDMwoMJlsH4m071R8TDRTJyRyQM2EU9Xu5UD1h43sp8hZprFl/0lSE2KC1umwY9Yn2G6bg37lp0be6A3H94HudeOoYbTSI5pxtpOQmUdcbPzgYrhgz14HibEzPA4D64uruO3pvX2tMq29ut8wPohERykuPyuH10tqebesIXIKAyF09ln5y4z/Y+1+938v2wpdzT7LNBw4yTW8GAw0t1u589kPWDAtja8sL/Brn8LYx6tOl2Oc82ZnYdPwxiE/vMyDPTCPiTMB9/UbTPD9lwvhlXudixC1NjKNnCWmiMkLToyBroHuuGJxDhUNHWwvazDuRl3NJhvvJVprDx7O/jlqQGR+p6ZnJPFed775pcZNgWDpS+axcE3fIptNs+XASc4ryiI2SDMV0zNT+In1alRHg3GG8gY3hYFgMs6W5DiS4kJj8Zdq9XB9GKmGVgEggXMAGMN4Z6eHkZoyuvysHGwafvLCd25TiQAAIABJREFUh0CEFAZC6NwT+vycPeic9z8HcSmmMt9PHHKN5nNvG3YwcO+GD2ls7+b+zy4cU9OzQnDxqV35GOX0nHTSk2J5zReds4ORGpjnnQvfeBMWfg62PgB//hi8+Rujq/zxROPlnu3ax8szpt326NY3u+Pj86eQGh/DkzsqBjhreK9zrmvtprWrx9lRw9ZrtuGHowZE5ndqekYSW1smmwJKd3KN0pdMa/EBXfner2yirrUrMDcNF/ItSezTMyjP+bQpYm8oG/5F1cXGiSptqtPiEbGiG4KORPf1Ye0elocTueMHQCinjPIzk5mZmcy+KmNL97W/74icAqNQuCdMzDNB+VE3Omdbr2m8UnRhQI1XHHKNZ23LhhwMvHmojnU7KvnqihmRI5kRIhKf25WPQaKjFCsKs3i9pNZM+fvCSA7ME9LgM7+HKx6D2gPw8h39+kow2kov9JWd1l6a2q1jTqoBkBgXzSdPn8qGfTW0ptgdF3zQOTscNQqyBhSyNx01BeV+Zpwj8Ts13ZJEZXs0vZbCwYFzV6upzyka7KYRHaX4aFF20I5jWnoiMVGKFyfdYNrZv+zFALO6eFBhIEBVYwc5IWxf/oD1qkE2g+06jgesV4XsGLwl9G12xhBT0xP7ZBquy4PN+uIqjjX2+zNWN3dy29N7AcbPTThvGRzabKZSB37Jj74F7XV+yzQcDJRrfOlrV7q9OXd093Lb03spyEzmptWFAe1PGB9cesa08fMd9UBKfDR1rV3MvH0DU9MTueXC2d5/JgvdfxeDxvxLjWe01cX/tseurxxm3ydbugBGdbvtobhycQ5PvHuMF47FcZWK8qkJSlmdKcB0kmoE4KjhINK+U3l2x5BTE08jvXqb8z2q7HUzUHDtFrj/JEvyJzIhyc921m6IiY5i2sREPjiVZLrgvnqvCdo9Fcx3tZrZlXnOzZ5tNuPhvGZ+8LLhw/FY61Iaorr5Qcw6pqp6qrWFB3qu5PmupdwdsqPwDsk4B8CXzskbtGykpozWbjqItTdyKonDQv4yEyDXurzn/c9DTALMWh3wLpzcNdzwq80lHGto577PnkZCbAS4mghChLO+uKqvuYEGqpo6RqQWJCBcCn/78EJf2efhPAYzzgCLctMpzE7hX8UnjLbch4zzkbo24qKjnJNJfYHz2Ek85GWYgUF18mxT0Drw/6lkE8SlwvRz+hZVNLRz4Pip4LhpuB6LJdm03T7nRkjLMYNCT83Dju813X5d9M21rV1099pG3lFjAFPTE3nOtpzl3cZmcHn3QzxnWx6R9SASOPuJ1prXS2qJj1FMTksY8SmjSKskDguOUfNAuYbNZgLnWashPnBf60HNUAawp7KJP207wueX5nL2jEEd4wVBcMPaTQfp7IlwS74AihCPj4GugUOhlOKKxTnsOtZEe0q+Txrnsto28ixJREcNmCGsK4Eki3urwVGKw6O6JMqu23bINbQ2Tblmnu/kIOKwlvWpzbaX5GUkUV7fho5NhNV3w/E98P4T7lf2WBhoZl9yQ6hxjkTtuickcPaT596v5q3D9dzxiXm8c/sqyu7/BG/eesGITR9FYiVxyJlYAKlTnQsEq3YaC6AAZRoOBjVDsWPttfHDf+8lMyWeWy+aG5R9CcJ4YFQM+gMoQjw5xjPOAJ85I4foKMWB7kyfLOnK6tqcCwPBZJwtYyfbDDAhMZb0pFiKrbnGn9oRkB7fa+5PLjKNzftPMjMrefBnEwTyLEmc6uyhqd1qvNBzlhjJkTvf8upi422e6hzAVzSEyMN5AJGoXfeEBM5+0NJp5acv7GdhzgS+8JHBco2RYDSNxkYMpYxc4+ib/Rfu/c+ZrmAuhReB4E6u8cjWI+yvaeEnly5gQmLwNGmCMNYZFYP+AIoQjzd3khAbRVri2C0ZykqN5/zZ2bxWl2Ys6TqG71fQa9McrW8fHBzWl44pmYaDvIwkDjfaIHtuf+Bcusk8DgicWzqtbC+rZ3UQ3TQG4vDMPtrQbu6ZF95n5CNv/Grwyo7CQBccGedQumqAc+fKkUxEBooEzn7wy5dKqGvt4qeXLnCeghpBRtNobETJW2YuAvWHTfC8/znTwjQxPWi7cJVrHK5t5cEtpVy0YDIXzo88axxBiGRGzaDfT3eg43YP57HeOfTKxTl80GmXqHlRIFjd1EF3r805cG5vgLbagAoDI5XplmSONrTB1EUmIHXINKYscsrobi2pxdqrR0TfDP2FikfrjaMJuUvgtCvg7Yeh6Vj/ip3NZhAzddGgbVQ2dpCVGi91PB6QwNlH9lU187e3y7n6I3kszAlesOYNo2U0NqIM1Dkf3wuN5UGTaTgYKNew2UxTm4SYKH58iXfNEARB6GfgoB+MPd1YGvSP1eYnrpw/J5umRO8t6fqs6AYGzvWHzOMYzThXN3XSO3kRtNdDzftQ+Z7bboETk2I5c/rEETmO3IwklMIUCDpYdZd53Hx3/7Ka982jm4xzRYg9nEcbEjj7gM2muWP9PjKS47j5YxGWLRkvWGYZs/byN0xRoIoybbaDzIzMZLaXNTDj9g28W9bAxadNITt17N8cBWEkcAz67/7UPHptmjOmhzbpMJKMxXbb7oiNjmLJokX0akX78eELO/s9nINrRRepTLck0WvTnEixN87ZutY4VhT2B849vTZeOXCS8+dkj9hsdUJsNJPTEih3ZJwB0nPh3G/Dvn/Dse1mmUNOMsWdVKMjpPrm0YYEzj7wrx0V7K5o4vaL5wbVe1HwAaVM1rn8TSPTyFsGyZlB3cX64ipeK3Hucvbs7qrIss8ShFHIyqIsALaV+tF+OwLRWnO8pZPJY9TD2ZXPLp1Jtc6k+siHw65bVtdGSnwMWSkDOirWlZialPTQ1AaFkjx7s5DmCtNfgQP/MYkdR5Yd2Hm0keYOK2tGSKbRdyyWJOeMM8Cy70DKZNh0m3Gjqi6G9OmQ7OwQ1WvTVDd1hNRRY7QhgbOX1Ld2cf/GAywtyOAzY2SKcdQSm2AqlWsPGLmGF929fGHtpoN0DbLPskWWfZYgjEIKMpOZlp7ItlI/2m9HGOuLqzj3/lfo7rHxr/eOjYuBddGkVOrjp9Fbdwg9jLPGEbujhpP2u/6QaTsdPfYKKfMsyXw66g0K372rf6G2wX++03eP2rz/BHHRUaywDyBH7FgykgcHzvEpsPou40S17ymPhYEnWjqx9mrJOA+BBM5e8rMXD9DW1cNPL10w5otAIpo962Df0/2/dzbB898OavA8KuyzBGEUopRiZVEmbx2qx9prG/4FEcr64ipue3ovNc3Geae5oyfymrqMEElTipjUU83equYh1yura3VjRVcyJvXNANmp8fwwdh0xNpfmWVZ7B0qMvvnsmRZS4kd24JCXmURdaxetXT3Of1j4OZgwHZ75mqkPOvLaoHtnZaPDik4yzp6QwNkLdpQ3sG5HJV9eUUDRpNRwH874Zss90OP5whQMRoV9liCMUlYUZnGqq4f3K5rCfSh+s3bTQTqszt3YIq6pywgxfdYC0lUbz7+zz+M6XT29VDZ2OAfOvVZTVDgG9c0AUVGKKare/R+bKzlc28qRujZWz80e8WNxWNIdc80673vKuFJp+6C1s3lQ4qmv+UmGZJw9IYHzMPT02rhj/T6mTkjg2xeMzZHyqMJTC1wvWuN6y6ixzxKEUciymZlEKdg6CnXOWmveLWugahzPSiVkm/vgh/t202l138r5WH07Wrs4ajQeBVvPmM04AzREe5BgTMgZ0W6BrkzPcLGkc7DlHujtcl7mknhyND+Zmj4+dPv+IIHzMDz6VjkHjp/izk/NJ3mEp1cELwigNa63iGe2IIwcE5JiOT03PeJ0zuuLq1h2/ysU3PoCy+5/xUl2UVbXxi9fOsjKta9y5f++jSex3riYlbLMBCCzu4pNdq97V9xa0dWVmMcxmnEGeD33m3ToOOeF9g6Umz88ydwpaX22jCOJw8u53DXj7EXiqbKxnUlp8cTHiIezJyQSHIKa5g5+9XIJ58/O4sL5Iz9KFLxg1Z1mask6ILPjZWtcX7j0jGkSKAvCCLGiMIuHXymlud0aEQ5FDs2yQ35R1dTBrU/vYVtpLYdr29hd0USUgmWzMvnu6iKsPTbufv5DJ7nGuJmVSs9Do1iYWM9TOyu5ZNHg66QjcM53FzhbZoXiKMNCa9Fn+GHJSX5leY7oU1UmobPqThpnXsqOJ17mxvND895TE2KxJMdxrMEl4zwhB5orBr9gQOKpsrGDXCkMHBIJnIfgp//ZT49N8+NPS0FgxODo5rXlHjNKtl+YvO3yJQhC+FlZmMlDW0p583AdF582JdyH41az3Gm18e9dVcyZnMrtF8/hkkXTnBqdxMdGs3bTQaqbOpiansgtF84eH4Pt2ATUhFzOjWvmp4fqqGxsH+TAUFbXRmZKHBMSBwyK6kuNB38Qu7xGGtMtSdxlW841l93CkvyMvuWv7qrEpkMj03CQZ0mivM4l4+xF4qmisZ3FeSPTnGWsIIGzB14vqeWFvTV8f00R0y0y+oooFl4pgbIgjGIW5aaTGh/DttLaiAicPWmTFfDiTSvd/m1cz0plFDCj/QRaw793VvGd1c66ZYcVnRN1pWNapgH9Xs5H69udAuct+0+SnRrPadMmhO5YLMm8W9bgvHCYxFNPr42a5k6xohuGgDTOSqkMpdTLSqlS++OgYYpSapFS6m2l1AdKqT1KqasC2Wco6LT2ctez+yjITOaG82aE+3AEQRDGFDHRUZw7y8LWkrph/YBDgTjp+IhlJvEt5SybZeGpXRXYbM7nsMxt4FwCmWNXpgGQMzGJKAXHBhTldffYeL2kllVzs4kaoW6B7sizJFHd3EFXj0sB58Ir4bv74O4m8zggCXW8pZNemxYrumEItDjwVmCL1roQ2GL/3ZV24Eta6/nAx4FfK6Uicq7GURwy5/+9SHl9OxfOnyQCeUEQhBFgRWEWVU0dHKlrG37lEearKwsGLRs3mmV/yJgBHY184bRUKho6eKes34btVKeV2lNdFGSm9K/fVg8djWM+4xwXE8WUCYkcbeiXSGwvq6e1q4fVIZRpgLGk07rfJcMbHOuKFd3QBBo4XwI8Zn/+GHCp6wpa6xKtdan9eTVwEhjZtjl+4CgOGWgz9NhbR8eFob0gCEKoWVlob79dEl53je4eG8/uriYuWpGdGi9OOt6QYZw1Vk9qIzUhhqd29LsyOHS1481Rw4Fru+st+0+SEBvFslmZIT0Oh8R0kCXdEDg8nCXjPDSBapwnaa1rALTWNUqpIZ29lVJLgTjgcID7DTpDGdrLxVMQBCG4TLckkW9JYltpHdctG5zxDRX/s2E/xcea+N0Xz4wIvfWoIMNIGONbyvn06XP5965K7r5kPmkJsRypawVgRtb4ctRwkGdJYtMHxrNZa83LH55g+awsEmJDO3vtaIIyyJJuCCobO1AKpkyQwHkohs04K6U2K6X2ufm5xJcdKaWmAH8Hrtdau+21qpS6QSm1Qym1o7Y2tFkIabMsCIIQWlYUZvH2kXq6e8LTfvvZ3VU8+lY5X1leIEGzL0zMBxQ0HOGKxbl0Wm28sKcGMBlnpfqbcADGUSM6HtKnh+VwQ8n0jGQa2ro51Wnl4IlTVDV1hKRboCsTk2JJjY9x0lsPR0VjO1PSEoiLkRYfQzHsp6O1Xq21XuDm51nghD0gdgTGJ91tQymVBrwA3KG1fmeIfT2itV6stV6clRVaNYcUhwiCIISWFYWZtHf3svNoY8j3XXLiFLf+ey9L8ifyw4vmhHz/o5rYBOPIUH+Y03MmUDQphXU7jD9wWV0rUyckOmdY60pNtjlq7NcM5Vn6nTU2f2gyzxeEIXBWSpGXmeRzxlkcNYYn0GHFc8C19ufXAs+6rqCUigOeAf6mtX4ywP2NGNJmWRAEIbScM9NCTJQKeRfBU51Wvv73nSTHx/DbL5xJbLRk2HwmowAajqCU4oqzcik+1sShk6coq2tzlmmA3Ypu7Ms0oD/Tfqyhnc37T3J6bjrZqeFpX52XkcyxBu8D56rGDtE3e0GgV4v7gTVKqVJgjf13lFKLlVJ/sq9zJbASuE4ptdv+syjA/QYdabMsCIIQWlITYjlz+kS2ldaFbJ9aa37w1B6ONrTz2y+cQXZaeIKaUU/GTGg4Apj7Z0yU4skdlYM9nHu6oLF8XBQGQn/GeefRRnZXNLEmDNnmgcdS0dBOT+/wUihrr42a5g5yxFFjWAIqDtRa1wOr3CzfAXzF/vz/gP8LZD+hYlwb2guCIISBFYWZ/HJzCfWtXVhS4kd8f3/aVsbGfcf50cVz+cgMy4jvb8ySMQM6GqCjkazUicyZnMojW4+gMdrxM6dPNPfThjLQveMmcE5NiCUjOY517xnpSii7BbqSb0mmx6apae4c1mKupqkTmxZHDW+Q+SlBEAQhbKwoykJreOPQyGedtx+p5/4XD3DRgsl8ZUX4nDzGBBZjSUfDEdYXV1FyohVHG5Tmjh5ue3qvsXOtL7WvPz6kGuuLqzjVaeVUVw/RSnGgpiVsx+KwpCv3okBQrOi8RwJnQRAEIWycNm0C6UmxIy7XONnSyY1PFJOXkcQDly9EqdB1cRuT2C3pqD/C2k0H6XaRAzjsXPs9nAsZ6zj6QVh7zRCiV2tuf2Zf2PpB+GJJV2EPnHOlOHBYJHAWBEEQwkZ0lGLZrEy2ldaOWPtta6+Nbz2+i9bOHv5wzVmkJsSOyH7GFRMLcFjSDWnnWlcKqVMgPjW0xxcGhuoHEQ6yU+OJj4nyypKusrGD6CjFlAmi+R8OCZwFQRCEsLKyMJMTLV2UnGgdke3/bOMB3itv5P7LTqNo0tgP4EJCbAKkTYOGw0PbudaVjotsM0ReP4ioKEWexTtLusrGDianJRAjDjPDIp+QIAiCEFaWO9pvj4At3Ya9NfzpjTKuPSePSxZJ8XdQscyAhiOe7Vw/VmQPnMdHYWAk9oOYnpHMMW+kGg3t5GaIvtkbJHAWBEEQwsq09ERmZiWzNUg65/XFVSy7/xUKbn2Bb/1jF3kZifzoE/OCsm1hABkzoP6wZzvXwljoagbL+Mg4R2I/iHxLEkcb2rDZhpZBSfMT7wnIjk4QBEEQgsHKoiwe336MTmuvc9c5H3EUaA3Ump5o6WLD3hqxGw02GTP7LOnc2rmWv2Eex4lUw/H+1246SHVTB1PTE7nlwtlh/b/Ly0ym02rj5KkuJnvQL3f19HLiVKc4aniJBM6CIAhC2FlZmMVf3yznvfIGVtilG/7grkCrs8fG2k0HJXAONg5njYYymDZx8N/7HDXGh1QDIq8fRF6GowV4m8fAubqpE63FUcNbRKohCIIghJ2PzMggNloFbEsXaQVaY5oBXs5uqTsEMYmmiFAICw5LuqND6JzFw9k3JHAWBEEQwk5SXAyL8zLYWhJYgWAkFmiNWSbmm0ePgXMJZM6CKAk1wsXU9ARiohRHGzxb0lU2mkGltNv2DvlvFgRBECKClUVZHDh+ipMtnX5v4/KzBmc3w12gNWaJtWeT6w+7/3tdybiSaUQiMdFR5ExMHNKSrqKhnZgoxeQ08XD2BgmcBUEQhIhgRWEmgN9yDa0120rrSI2PZsqEBGeHhwjSnY4pMma4zzhbO6Hp2Lhx1Ihk8ixDW9JVNppCxugo6abpDVIcKAiCIEQE86akYUmOY1tpLZedlePz61/cd5xdx5r42WWncdWS6SNwhMIgMmbAgf8MXt5wGNDjxlEjksmzJLHrWCNaa7et5isb20Xf7AOScRYEQRAigqgoxfLCTN44VDes76wr3T027n/xALMnpXL5WbkjdITCICwzob0eOpqcl/c5akjgHG7yLMmc6uyhsd3q9u8VjR3iqOEDEjgLgiAIEcPKwizqWrv5sKbFp9f9Y/tRjta3c+vFc2TKOZT0WdK5yDXqDplHy6zQHo8wiIGWdK50WnupPdUlGWcfkMBZEARBiBj80Tk3d1h5cEspy2dl8tEi/z2gBT/I8GBJV1cCE3IhLjn0xyQ4kZ/pCJwH65yrmhyOGhI4e4sEzoIgCELEkJ2WwJzJqWwr9d6W7nevHaK5w8ptF89xq+EURhBPlnR1JZJtjhByJiahlPvA2WFFJ1IN75HAWRAEQYgoVhRmsqO8kfbunmHXrWxs569vlvPZM3KYP3VCCI5OcCIuyVjSDQyctYb6Q2JFFyEkxEYzJS3BrVSjosHR/EQCZ2+RwFkQBEGIKFYWZdHda2P7kYZh1/35poMo4OYLJUgLG66WdKdqoLtVCgMjiDxLMuVuAufKxg7ioqPITo0Pw1GNTiRwFgRBECKKJfkZxMdEsXUYucbeymbW767mKysKmDJBNJphI2OGcxOUulLzKIFzxJBnSeJYgzupRjvTJiYSJQW1XiOBsyAIghBRJMRGs7QgY8gCQa019274EEtyHF8/b2YIj04YRMYMaK+Dzmbze58VncwCRAp5lmTqWrtp7XKWP1U0doijho8EFDgrpTKUUi8rpUrtjxOHWDdNKVWllHo4kH0KgiAIY5+VhVkcOtlKtb3q35VXDpzknSMN3LS6kNSE2BAfneCEqyVdXSnEpUDqlPAdk+BEnsW9JV2VND/xmUAzzrcCW7TWhcAW+++e+AnweoD7EwRBEMYBK+22cu7cNXp6bfzPhv3MyEzmc0ulQ2DYsbhY0tWXGkcNcTiJGPoD5365Rkd3L3Wt3VIY6COBBs6XAI/Znz8GXOpuJaXUWcAk4KUA9ycIgiCMA4ompZCdGs9WN3KNf+2o4HBtGz+8aA6x0aI4DDsTC8xj/YCMs8g0Ioo8i/HTHhg4VzY6HDUk4+wLgV5xJmmtawDsj9muKyilooBfALcEuC9BEARhnKCUYkVhFm8eqqN3QPvt1q4efvVyCUvzM/jYvElhPEKhj7gkSJ1qMs7dbdBcIYWBEUZKfAyZKXFOUg2Hh7NknH1j2MBZKbVZKbXPzc8lXu7jm8AGrXWFF/u6QSm1Qym1o7bWe/N7QRAEYeyxsiiTpnYr+6qa+5Y98vph6lq7uf0Tc6XZSSSRMQMaDve7a0jgHHG4WtI5Ms650jXQJ2KGW0FrvdrT35RSJ5RSU7TWNUqpKcBJN6udA6xQSn0TSAHilFKtWutBemit9SPAIwCLFy/Wrn8XBEEQxg/LZ5n221tLajk9N53jzZ08su0Inzp9Koty08N8dIITlhlwcKM4akQweRlJvHOkvu/3isYO4mOiyEoRD2dfCFSq8Rxwrf35tcCzritorb+otZ6utc4Hbgb+5i5oFgRBEISBWFLiWTAtrc+W7pcvH8Rmgx9cODvMRyYMImMGtNVCdTGg+p02hIghz5JMTUsnndZeoN/DWWZufCPQwPl+YI1SqhRYY/8dpdRipdSfAj04QRAEYXwzKS2Bd8sbKLj1BdbtqGTZrAxyM0STGXFk2J01Sl6E9OkQK9P/kUaeJQmt+yUalY0d5Iq+2WcCCpy11vVa61Va60L7Y4N9+Q6t9VfcrP+o1vrGQPYpCIIgjA/WF1f1ZZsd2r23jzSwvrgqfAcluMeRYa4/JDKNCMVhSVdeZwLnigbxcPYH8fERBEEQIpK1mw7S3WNzWtZptbF208EwHZHgkYyC/udSGBiROCzpyuvbaO3qobHdKo4afiCBsyAIghCReOoa6Gm5EEYOvADKHlK8/wTsWRfe4xEGMTEpltSEGI41tFNlt6ITRw3fkcBZEARBiEimpru/qXtaLoSJPevg+W+Dts8OdDSa3yV4jiiUUuRbkimvb6eiwdH8RDLOviKBsyAIghCR3HLhbBJjo52WJcZGc4u4akQWW+4Bq8ssgLXDLBciiumWJI7Vt0nXwACQwFkQBEGISC49Yxr3ffY0pqUnooBp6Ync99nTuPSMaeE+NGEgzZW+LRfCRr4licrGDsrr20mMjcaSHBfuQxp1DNsARRAEQRDCxaVnTJNAOdKZkGPabLtbLkQUeRnJ9Ng028sayBEPZ7+QjLMgCIIgCP6z6s7Bvs2xiWa5EFE4LOn217SITMNPJHAWBEEQBMF/Fl4Jn3oIJuQCyjx+6iGzXIgo8jOT+55LYaB/iFRDEARBEITAWHilBMqjgOzUeBJio+i02sSKzk8k4ywIgiAIgjAOeHZ3NT29pg/nH147LF04/UACZ0EQBEEQhDHO+uIqbnt6Lz02Ezg3tFu57em9Ejz7iATOgiAIgiAIY5y1mw7SYe11WtZh7ZUW9j4igbMgCIIgCMIYR1rYBwcJnAVBEARBEMY40sI+OEjgLAiCIAiCMMaRFvbBQezoBEEQBEEQxjiODpxrNx2kuqmDqemJ3HLhbOnM6SMSOAuCIAiCIIwDpIV94IhUQxAEQRAEQRC8QAJnQRAEQRAEQfACCZwFQRAEQRAEwQskcBYEQRAEQRAEL5DAWRAEQRAEQRC8QGmtw30MblFK1QJHw30cEUomUBfugxCGRc7T6EDO0+hAztPoQM7T6EDO02DytNZZw60UsYGz4Bml1A6t9eJwH4cwNHKeRgdynkYHcp5GB3KeRgdynvxHpBqCIAiCIAiC4AUSOAuCIAiCIAiCF0jgPDp5JNwHIHiFnKfRgZyn0YGcp9GBnKfRgZwnPxGNsyAIgiAIgiB4gWScBUEQBEEQBMELJHCOAJRSuUqpV5VS+5VSHyilvmNfnqGUelkpVWp/nGhfrpRSDymlDiml9iilzrQvX6SUetu+jT1KqavC+b7GGsE6TwO2l6aUqlJKPRyO9zNWCeZ5UkpNV0q9ZN/Wh0qp/PC8q7FHkM/TA/Zt7Levo8L1vsYifpyrOfZ7UZdS6maXbX1cKXXQfh5vDcf7GYsE6xx52o4wAK21/IT5B5gCnGl/ngqUAPOAB4Bb7ctvBX5mf34xsBFQwNnAdvvyIqDQ/nwqUAOkh/v9jZWfYJ2nAdt7EHgceDjc720s/QTzPAGvAWvsz1OApHC/v7HyE8Tr3rkitsOxAAADuUlEQVTAm0C0/edt4KPhfn9j6cePc5UNLAHuBW4esJ1o4DAwA4gD3gfmhfv9jYWfIJ4jt9sJ9/uLpB/JOEcAWusarfUu+/NTwH5gGnAJ8Jh9tceAS+3PLwH+pg3vAOlKqSla6xKtdal9O9XASWBYM2/BO4J1ngCUUmcBk4CXQvgWxgXBOk9KqXlAjNb6Zfu2WrXW7aF8L2OZIH6fNJCACcTigVjgRMjeyDjA13OltT6ptX4PsLpsailwSGt9RGvdDfzTvg0hQIJ1jobYjmBHAucIwz4VfAawHZikta4B88+MGSGC+SeuGPCySlz+sZVSSzE3ksMje8Tjk0DOk1IqCvgFcEuojne8EuD3qQhoUko9rZQqVkqtVUpFh+rYxxOBnCet9dvAq5gZthpgk9Z6f2iOfPzh5bnyxLD3LiFwAjxHnrYj2JHAOYJQSqUA/wZu0lq3DLWqm2V99ij2LMzfgeu11rbgHqUQhPP0TWCD1rrCzd+FIBGE8xQDrABuxkxpzgCuC/JhjnsCPU9KqVnAXCAHE4RdoJRaGfwjFXw4Vx434WaZWHsFkSCco6BuZywigXOEoJSKxfyT/kNr/bR98YkBU/tTMNILMKP03AEvzwGq7eulAS8Ad9inM4UgEqTzdA5wo1KqHPg58CWl1P0hOPxxQ5DOUyVQbJ9W7gHWA04FnkJgBOk8fQZ4xy6lacXooM8OxfGPJ3w8V57weO8SAidI58jTdgQ7EjhHAPYK8D8D+7XWvxzwp+eAa+3PrwWeHbD8S/Yq87OBZq11jVIqDngGowN8MkSHP24I1nnSWn9Raz1da52PyWb+TWst1eVBIljnCXgPmKiUctQJXAB8OOJvYJwQxPN0DDhPKRVjv+Gfh9FlCkHCj3PlifeAQqVUgf1+9Tn7NoQACdY5GmI7gh1pgBIBKKWWA9uAvYBDWnE7Rle0DpiOuTlcobVusP9jPwx8HGjHSDJ2KKWuBv4KfDBg89dprXeH5p2MbYJ1nly2eR2wWGt9Y0jexDggmOdJKbUGo0dXwE7gBntRkxAgQbzuRQO/A1Zipv1f1Fp/L6RvZozjx7maDOwA0uzrt2KcGVqUUhcDv8Y4bPxFa31vSN/MGCVY5whY6G47WusNIXorEY8EzoIgCIIgCILgBSLVEARBEARBEAQvkMBZEARBEARBELxAAmdBEARBEARB8AIJnAVBEARBEATBCyRwFgRBEARBEAQvkMBZEARBEARBELxAAmdBEARBEARB8AIJnAVBEARBEATBC/4/Q0Qx7+bAx2MAAAAASUVORK5CYII=\n",
      "text/plain": [
       "<Figure size 864x288 with 1 Axes>"
      ]
     },
     "metadata": {
      "needs_background": "light"
     },
     "output_type": "display_data"
    }
   ],
   "source": [
    "#com+ your code here \n",
    "ar_result = ar.fit(maxlag=optlag)\n",
    "prediction = ar_result.predict(start=optlag)\n",
    "\n",
    "plt.figure(figsize=(12, 4))\n",
    "plt.plot(time_series, '-o', label='true')\n",
    "plt.plot(prediction, '-o', label='model')\n",
    "plt.legend(fontsize=12);\n"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "As for the linear regression, we may want to look at the \"size\" of the residuals. For example you could display the Mean Absolute Error (MAE) using `mean_absolute_error` and feeding the original time series from `optlag` onwards and compare it with the prediction."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "#com+ add your code here to show the MAE\n",
    "print('MAE = {0:.3f}'.format(mean_absolute_error(time_series[optlag:], prediction)))\n"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "### Autoregression with Sklearn\n",
    "\n",
    "Autoregression can also be implemented using `sklearn`. However, it doesn't provide direct support to handle time series, which means that you have to reorganise the data before estimating the model:"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "collapsed": true
   },
   "outputs": [],
   "source": [
    "def organize_data(to_forecast, window, horizon=1):\n",
    "    \"\"\"\n",
    "     Input:\n",
    "      to_forecast, univariate time series organised as numpy array\n",
    "      window, number of observations to use in a forecast window\n",
    "      horizon, horizon of the forecast\n",
    "     Output:\n",
    "      X, a matrix where each row contains a forecast window\n",
    "      y, the target values for each row of X\n",
    "    \"\"\"\n",
    "    shape = to_forecast.shape[:-1] + (to_forecast.shape[-1] - window + 1, window)\n",
    "    strides = to_forecast.strides + (to_forecast.strides[-1],)\n",
    "    X = np.lib.stride_tricks.as_strided(to_forecast,\n",
    "                                        shape=shape,\n",
    "                                        strides=strides)\n",
    "    y = np.array([X[i+horizon][-1] for i in range(len(X)-horizon)])\n",
    "    return X[:-horizon], y\n",
    "\n",
    "\n",
    "X, y = organize_data(time_series, optlag)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "We can now fit an autoregressive model which simply corresponds to a `LinearRegression` now:"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "#com+ add your code here to fit a linear regression and show the results\n",
    "lr = LinearRegression()\n",
    "lr.fit(X, y)\n",
    "prediction = lr.predict(X)\n",
    "\n",
    "plt.figure(figsize=(12, 4))\n",
    "plt.plot(time_series.values, '-o',)\n",
    "plt.plot(np.arange(optlag, len(time_series)), prediction, '-o', label='prediction')\n",
    "plt.legend(fontsize=12);\n"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "The results are the same, yay!\n",
    "\n",
    "***TIP***: Now that you know how to implement autoregression with sklearn, you're also able to create custom autoregressive models using other regressors implemented in sklearn."
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## ARMA model\n",
    "\n",
    "AR models a point in the time series as a linear model of the previous values. The mismatch $e_t$ is assumed to be \"noise\".\n",
    "However there could still be information in the series of $e_t$! How about we add the past errors as additionnal features?\n",
    "\n",
    "This leads to the **ARMA** model with an autoregressive part that you will recognise and a part that corresponds to a moving average:\n",
    "\n",
    "$$\n",
    "y_{t} = c + \\underbrace{ \\phi_{1}y_{t-1} + \\phi_{2}y_{t-2} + \\dots + \\phi_{p}y_{t-p} }_{AR(p)} + \\underbrace{ \\theta_{1}e_{t-1} + \\theta_{2}e_{t-2} + \\dots + \\theta_{q}e_{t-q} }_{MA(q)} +e_{t},\n",
    "$$\n",
    "\n",
    "ARMA models are also implemented in `statsmodels` and their implementation is consistent with the one of AR models. \n",
    "\n",
    "* Create an ARMA model with `tsa.ARMA`, specify $p=3$ and $q=3$\n",
    "* Fit and predict\n",
    "* Display\n",
    "\n",
    "Since the result will look almost identical to just using AR, you will want to show the MAE as well. "
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "#com+ add your code here to fit an ARMA model\n",
    "arma = tsa.ARMA(time_series, order=(3, 3))\n",
    "arma_result = arma.fit()\n",
    "prediction = arma_result.predict(start=3)\n",
    "\n",
    "plt.figure(figsize=(12, 4))\n",
    "plt.plot(time_series, '-o', label='true')\n",
    "plt.plot(prediction, '-o', label='model')\n",
    "plt.legend();\n"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "#com+ add your code here to show the MAE\n",
    "print('MAE = {0:.3f}'.format(mean_absolute_error(time_series[3:], prediction)))\n"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## Forecasting\n",
    "\n",
    "Until now we have tuned the model to *fit* the time series. \n",
    "This makes us able to find a description of the data, but doesn't give forecast.\n",
    "\n",
    "**Note**: forecasting with time series data is **tricky** and usually basic methods do not really provide very good results (especially on realistic data). ARMA models are nice because they are simple but do not expect fantastic performances. (On the other hand, predicting the future is hard! -- who would have thought).\n",
    "\n",
    "* Separate the time series into a training set and a test set formed of the last 8 points. \n",
    "* Fit an AR model on the training data and try to find the optimal lag using the `BIC` criterion\n",
    "* fit the model, predict and show the prediction on the original time series. Did it do a good job? \n",
    "* compute the MAE"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "#com+ your code here\n",
    "train = time_series[:-8]\n",
    "test = time_series[-8:]\n",
    "\n",
    "ar = tsa.AR(train)\n",
    "optlag = ar.select_order(10, ic='bic', method='mle')\n",
    "arfit = ar.fit(maxlag=optlag)\n",
    "\n",
    "prediction = arfit.predict(end=len(time_series))[-len(test):]\n",
    "\n",
    "plt.figure(figsize=(8, 6))\n",
    "plt.plot(time_series.values, '-o', label=\"original data\")\n",
    "plt.plot(prediction, '--o', label='prediction')\n",
    "\n",
    "plt.legend(fontsize=12)\n",
    "\n",
    "print('Out of sample MAE = {0:.3f}'.format(mean_absolute_error(test, prediction)))\n"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "<h2>Cross-Validation for Time Series Forecasting (demo) </h2>\n",
    "\n",
    "In this procedure, there is a series of test sets, each consisting of a single observation. The corresponding training set consists only of observations that occurred prior to the observation that forms the test set. Thus, no future observations can be used in constructing the forecast. \n",
    "\n",
    "Let's build a diagram to illustrate the procedure:"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "test_size = 4\n",
    "\n",
    "plt.figure(figsize=(15, 10))\n",
    "for i, split in enumerate(range(len(time_series)-test_size, len(time_series))):\n",
    "    plt.subplot(test_size, 1, i+1)\n",
    "    plt.plot(time_series.values[:-test_size+i], \n",
    "             '-o', label='train set at fold %d' % (i+1))\n",
    "    plt.plot(range(len(time_series)-test_size-1+i, len(time_series)), \n",
    "             time_series.values[-test_size-1+i:], \n",
    "             '--', color='gray', label='remaining data')\n",
    "    plt.plot(split, time_series.values[split], \n",
    "             'o', label='test point at fold %d' % (i+1))\n",
    "    plt.xlim([0, len(time_series)])\n",
    "    plt.legend()"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "We can now implement the cross-validation training and testing the autoregressive model at each step:"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "collapsed": true
   },
   "outputs": [],
   "source": [
    "start_split = 4*8\n",
    "predictions = []\n",
    "\n",
    "for split in range(start_split, len(time_series)):\n",
    "    train = time_series[:split]\n",
    "    test = time_series[split:]\n",
    "    ar = tsa.AR(train)\n",
    "    ar_result = ar.fit(maxlag=4)\n",
    "    prediction = ar_result.predict(end=len(time_series))[-len(test):]\n",
    "    predictions.append(prediction.values[0])"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "len(predictions), len(time_series)-start_split\n",
    "plt.plot(time_series.values, '-o', label='true')\n",
    "plt.plot(range(start_split, len(time_series)), predictions, \n",
    "         '-o', label='out of sample prediction')\n",
    "plt.legend();"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "print(mean_absolute_error(time_series.values[start_split:], predictions))"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Wrap up\n",
    "- Using lags of the target variable as feature one can build an autoregressive model\n",
    "- Adding moving average terms and differencing lead to ARIMA\n",
    "- The order of the ARIMA models can be hard to tune\n",
    "- Time series models can be used to describe the data and make forecast\n",
    "- Cross validation can be used to estimate the error of the forecast"
   ]
  }
 ],
 "metadata": {
  "kernelspec": {
   "display_name": "Python 3",
   "language": "python",
   "name": "python3"
  },
  "language_info": {
   "codemirror_mode": {
    "name": "ipython",
    "version": 3
   },
   "file_extension": ".py",
   "mimetype": "text/x-python",
   "name": "python",
   "nbconvert_exporter": "python",
   "pygments_lexer": "ipython3",
   "version": "3.7.0"
  }
 },
 "nbformat": 4,
 "nbformat_minor": 2
}
